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Abstract 

Noncooperative game-theoretic tools have been increasingly used to study many important resource 
allocation problems in communications, networking, smart grids, and portfolio optimization. In this paper, 
we consider a general class of convex Nash Equilibrium Problems (NEPs), where each player aims to solve 
an arbitrary smooth convex optimization problem. Differently from most of current works, we do not 
assume any specific structure for the players' problems, and we allow the optimization variables of the 
players to be matrices in the complex domain. Our main contribution is the design of a novel class of 
distributed (asynchronous) best-response- algorithms suitable for solving the proposed NEPs, even in the 
presence of multiple solutions. The new methods, whose convergence analysis is based on Variational 
Inequality (VI) techniques, can select, among all the equilibria of a game, those that optimize a given 
performance criterion, at the cost of limited signaling among the players. This is a major departure from 
existing best-response algorithms, whose convergence conditions imply the uniqueness of the NE. Some of 
our results hinge on the use of VI problems directly in the complex domain; the study of these new kind of 
Vis also represents a noteworthy innovative contribution. We then apply the developed methods to solve 
some new generalizations of SISO and MIMO games in cognitive radios and femtocell systems, showing a 
considerable performance improvement over classical pure noncooperative schemes. 



X 

^ 1 Introduction and Motivation 

In recent years, there has been a growing interest in the use of noncooperative games to model and solve 
resource allocation problems in communications and networking, wherein the interaction among several agents 
is by no means negligible and centralized approaches are not suitable. Examples are power control and 
resource sharing in wireless/wired peer-to-peer networks, cognitive radio systems (e.g., [UElISlSlISllllITlElEl 
Uni [m [E]), distributed routing, flow and congestion control, and load balancing in communication networks 
(e.g., |13 | [T H [T5] and references therein), and smart grids (see |16| 117) and references therein). Two recent 
special issues on the subject are |18| 119) . 
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Among the variety of models and solution concepts proposed in the literature, the Nash Equilibrium Prob- 
lem (NEP) plays a central role and has been used most to model interactions among individuals competing 
selfishly for scarce resources. In a NEP there is a finite number / of players; each player i makes decisions 
on a set of variables Xj belonging to a given feasible set Xj S Qi- The goal of each player i is to minimize his 
own objective function fi{'Ki,:x.-i) over Qi while anticipating the reactions x_j = (xj)j_^j^-^ from the rivals: 

minimize fii^i, x_i) 

(1) 

subject to Xj G Qj. 

The NEP is the problem of finding a vector x* = (x*)^^-^ such that each x* belongs to Qi and solves the 
player's problem (given x* J: 

/,(x^ X* J < /,(x,, x:^J, Vx, GQj. (2) 

Such a point x* is called a Nash Equilibrium (NE) or, more simply, a solution of the NEP. In words, a NE 
is a feasible strategy profile x* such that no single player can benefit from a unilateral deviation from x*. 

In this paper we focus on NEPs in the general form ([1]), in the following setting: i) the optimization 
variables of each player can be either real vectors or complex matrices; ii) each optimization problem in 
([1]) is convex for any given feasible x_j; and iii) players' objective functions are continuously differentiable 
in all the variables (more precisely, functions of complex variables are assumed to be M-differentiable, see 
Sec. [6]). We will term such a game (real or complex) player-convex NEP. Note that assumptions ii) and 
iii) are mild and quite standard in the literature, see for example [Tt[2|l3lHl[5|[6l[71[8l[9l[T0l[TT | [T8 | [T9] 
where special instances of the player-convex NEP ([1]) are studied. The convexity assumption ii) makes the 
NEP numerically tractable (a NE may not even exist otherwise) while, to date, the differentiability of the 
players' functions seems indispensable to analyze distributed solution methods |20| [2T] (unless the game has 
a very specific structure, like in potential or supermodular games; see, e.g., \22\ I23| and references therein). 
Motivated by recent applications of noncooperative models in MIMO communications [71 [HI [TOl [TT| [24], we 
also allow, according to i), players' optimization variables to be complex matrices, which significantly enlarges 
the range of applicability of model ([T|). To the best of our knowledge, this is the first work where a complex 
NEP in the general form ([1]) is considered. 

While the solution analysis (existence, etc.) of a real player-convex NEP relies on standard results in 
game theory (see, e.g., the seminal work [25], or |20) for more recent results), the development of distributed 
solution algorithms is much more involved. The goal of this paper is to address this difficult task in the broad 
setting described above. We are interested in the design and analysis of (possibly) asynchronous iterative 
best-response algorithms, suitable for solving real and complex player-convex NEPs, even in the presence of 
multiple NEs. By "best-response" algorithms we mean iterative schemes where the players iteratively choose 
the (feasible) strategy that minimizes their cost functions, given the actions of the other players; the reason 
for our emphasis on best-response schemes will be described shorty. 

1.1 Literature review 

The study of iterative algorithms for (special cases of) player-convex NEPs has been addressed in a number of 
papers, under different settings and assumptions; the main features and limitations of current state-of-the-art 
approaches are discussed next. 
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A first class of papers is composed of works motivated by specific applications, some examples are [H 
[21 131 m El [71 m [9l [TOl where different resom'ce allocation problems in communications are modelled as 
noncooperative games and solved via iterative algorithms; all these formulations are special cases of the NEP 
([T|). A key feature of all these models is that the best-response of each player (i.e., the optimal solution of each 
player's optimization problem) is unique and can be expressed in closed form; this simplifies enormously the 
application of standard fixed-point arguments to the study of the convergence of best-response algorithms. 
The same line of analysis is developed also in |26| and |27) for quadratic and non-quadratic (but with only 
two-players) real convex NEPs, respectively. A monotonicity-based approach is instead used in |28|I29|. Even 
though algorithms in [261 [27l [28] do not require a closed form solution of players' optimization problems, 
they can be computationally very demanding and the convergence conditions are based on assumption whose 
verification for games arising from realistic applications remains elusive. Last but not least, convergence 
conditions of the algorithms proposed in all the aforementioned papers imply the uniqueness of the NE. 

A more general and powerful methodology suitable for studying noncooperative games is offered by the 
theory of finite-dimensional Variational Inequalities (Vis) |30| . VI and complementarity problems have a long 
history and have been well documented in the literature of operation research |30| . but only recently they have 
been brought also to the attention of the signal processing, communications, and networking communities 
[21 H [g [ini [311 [32|. Given a subset K of M" and a vector- valued function F : /C M", the VI problem, 
denoted by VI (/C, F), consists in finding a point x* G /C such that 

(x - x*)^ F(x*) > Vx G /C. (3) 

The VI approach to noncooperative games hinges on an easy equivalence between real player-convex NEPs 
as in ([H) and the (partitioned) VI problem VI(/C,F) in ([3]), with K, = 0^=1 Qi and F = (Vxji(x))f^i 
(intended to be a column vector), where Vxj/i(x) denotes the gradient of fi with respect to x^. Based 
on this equivalence, one can solve a real player-convex NEP by focusing on the associated VI problem and 
taking advantage of the many (centralized and distributed) solution methods available in the literature for 
partitioned Vis ^ Vol. II]. 

In the effort of obtaining distributed schemes for NEPs, researchers have focused on so called projection 
algorithms for partitioned Vis [301 Ch. 12]. When applied to NEPs, projection algorithms work as follows: 
at each iteration, all the players choose simultaneously their own strategies employing a (projected) gradient 
update. This approach has been followed in a number of works, e.g., |33| [Ml [35]; continuous-time strategy 
evolutions based on players' gradient responses have been studied also in |25| [36] (without using the VI 
formalism). Even though algorithms in |33[ [M] [35] seem to be a most natural choice to solve NEPs while 
using existing VI solution methods, they suffer, together with |25| 136). from some drawbacks, which strongly 
limit their applicability in practice. First, they are not "incentive compatible", meaning that selfish users 
may deviate from them, unless they are imposed by some authority as a protocol to follow. Second, and 
most importantly, they generally converge very slowly; this has been observed in a number of different 
applications (see, e.g., numerical results in |34[ I35| [55]). For instance, our experiments show that gradient- 
response algorithms need about two orders of magnitude more iterations than best-response schemes to reach 
comparable performance (see. Figure [6] in Sec. [8]); this makes them not appealing in practice, especially for 
the design of wireless systems, where iterations correspond to physical transmissions. 
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A different approach to the design of algorithms for partitioned Vis has been followed in [371 138] , where 
the authors investigated the local and global convergence of various iterative synchronous methods that 
decompose the original VI problem into a sequence of simpler lower-dimensional VI subproblems. When 
the VI problem comes from a real player-convex NEP, some of the iterative algorithms in \37\ [38] can be 
interpreted as synchronous best-response schemes solving the NEP. Unfortunately the convergence analysis 
in |37| I38| . based on contraction arguments, leads to abstract convergence conditions, whose verification 
in practice seems not possible. Easier conditions to be checked have been obtained recently in |20) for 
simultaneous best-response algorithms, still using the VI approach. However, conditions in |201 l37l [38] are 
applicable only to a restricted class of real NEPs; they indeed imply the (uniformly) strongly convexity of the 
players' cost functions and the uniqueness of the NE. In the presence of multiple solutions, the distributed 
computation of even a single NE of real/ complex NEPs via best-response algorithms becomes a difficult and 
unsolved task. 

The analysis of the current literature carried out so far leads to the following conclusions: When it comes 
to distributed computation of NE via best-response dynamics, the following issues arise: i) the convergence 
analysis and algorithms apply only to a restricted class of NEPs, whose players' cost functions and feasible 
sets have a very specific structure, leaving outside player-convex NEPs in the general form ([1]); ii) the best- 
response mapping of each player must be unique and/or is required to be computed in closed form; iii) 
convergence is obtained only under conditions implying the uniqueness of the NE; and iv) none of current 
results and Vl-based methodologies can be applied to study and solve complex player-convex NEPs, which 
arise naturally, e.g., from applications in MIMO communications. 

1.2 Main contributions 

In order to address the key issues listed at the end of the previous subsection, in this paper we introduce 
several new developments that are summarized next. 

1. Building on our recent contributions |20 | l2 H [3 H l39). we develop a Vl-based unified theory for the study 
and design of distributed best-response algorithms for the solution of real player-convex NEPs, having 
(possibly) multiple solutions. Our unified framework has many desirable properties, such as: 

— It provides a systematic methodology for analyzing old and new algorithms, simplifying greatly the 
application of game-theoretical models to new problems; 

— It includes (both synchronous and) totally asynchronous methods, thus improving on the more 
traditional synchronous methods studied in the literature, see e.g. [20l [25l [33l (Ml (Ml [36])- Note that 
in many applications asynchronous algorithms play a fundamental role and are the only methods that 
can reasonably be used; 

— It advances the state-of-the-art algorithms by relaxing the convergence conditions. In particular, the 
proposed distributed best-response schemes are suitable to solve player-convex NEPs having multiple 
solutions; moreover, if some (limited) signaling among the players is allowed, one can also control the 
quality of the achievable solution by forcing convergence to a NE that optimizes a further performance 
criterion (i.e., performing an equilibrium selection). This feature is very appealing in the design of 
practical wireless systems, where algorithms with unpredictable performance (like [Il[2l[3llU[5l[71l8l[9l 



4 



[inilll] and \25 \ l33 l [M l I35 | 136), when multiple solutions are present) are not acceptable; 

— It does not require the players' best-response to be unique or given in closed form; 

— It allows us to clearly gauge the trade-off between signaling and characteristic of the resulting 
algorithms. 

2. We develop an entirely new theory for the study of Vis in the complex domain, which allows us to 
extend all the results we established in the real case to player-convex NEPs whose players' optimization 
variables are complex matrices. Key to this extension is the development of several new technical results 
that are of independent interest and useful in many other applications. Among these we single out 

— A Taylor expansion of real-valued functions of complex matrices that is amenable to our MIMO 
applications; 

— The minimum principle for constrained convex optimization problems in the domain of complex 
matrices, generalizing the already known complex gradient-vanishing conditions obtained in |40| for the 
unconstrained minimization of real- valued functions of complex variables; and 

— A new definition of Hessian matrix for real-valued functions of complex matrices (and Jacobian 
matrix of complex-valued VI functions) along with a relaxed concept of positive (semi-)definiteness, 
which are instrumental to establish the connection with the convexity properties of the function (and 
the monotonicity/P properties of the complex- valued VI function). 

To the best of our knowledge the above features constitute a substantial advancement in the distributed so- 
lution methods of noncooperative games, which enlarges considerably scope and flexibility of game-theoretical 
models in wireless distributed (MIMO) networks. We believe that the results in this paper bring the use of 
game-theoretic methods one step closer to "technological" maturity and will facilitate further developments 
and new applications in several fields. 

In order to illustrate our techniques we consider some new MIMO games over vector Gaussian Interference 
Channels (ICs), modeling some distributed resource allocation problems in SISO and MIMO CR systems and 
femtocells. These games are examples of NEPs that cannot be handled by current methodologies. The 
application of our framework leads to a unified algorithm converging to different types of solutions according 
to the selected degree of cooperation and signaling (in the form of pricing at the level of players' objective 
functions). The choice of the specific (pricing) scheme will depend then on the trade-off between signaling 
and performance that the users are willing to exchange/ achieve. Numerical results show the superiority of 
our approach with respect to plain noncooperative solutions as well as good performance with respect to 
centralized solutions, in favor of very limited signaling among the players. 

The paper is organized as follows. Sec. [2]introduces the just mentioned new resource allocations problems. 
Sec. [3] presents some results on partitioned Vis geared towards the analysis of NEPs; special emphasis is 
given to some classes of vector functions F and its properties that play a key role in the solution analysis 
of the VI as well as convergence analysis of distributed algorithms for NEPs. Establishing the connection 
between partitioned Vis and real convex-player NEPs and building on results introduced in Sec. [3l Sec. S] 
focuses on the solution analysis of real convex-player NEPs (more general results are given in Appendix [A]) . 
Sec. [5] and Sec. [6] constitute the core theoretical part of the paper; in Sec. [5] we provide various distributed 
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algorithms for solving real player-convex NEPs in several significant settings along with their convergence 
properties; Sec. [6] generalizes the main results in Sec. [3l HI and [5] on real Vis to the complex case. Sec. 
[7] shows how to apply the developed machinery to the concrete resource allocation problems introduced in 
Sec. [21 whereas Sec. [8] provides some numerical results corroborating our theoretical findings. Finally, Sec. 
[9] draws some conclusions. The proof of the major results of the paper is given in Appendix lAlHl 

2 Motivating Examples: Noncooperative Games Over Gaussian ICs 

To motivate and illustrate our new results more in detail, we start introducing some novel resource allocation 
problems over SISO frequency-selective and MIMO Gaussian ICs, widely extending formulations that have 
already been studied in the literature. We will show that these problems cannot be analyzed and solved using 
current results and algorithms, but call for a more general theory. 

The IC is suitable to model many practical multiuser systems, such as digital subscriber lines, wireless 
ad- hoc and Cognitive Radio (CR) networks, peer-to-peer systems, multicell OFDM/TDMA cellular systems, 
and Femtocell-based networks. To have specific case studies at hand, we will focus on CR systems and 
Femtocell-based networks; the proposed techniques can be readily applied also to the other aforementioned 
network models. 



We consider an /-user A^-parallel Gaussian interference channel, modeling a CR system composed of Q 
secondary users (SUs) and P primary users (PUs). In this model, there are / transmitter-receiver pairs— the 
SUs— where each transmitter wants to communicate with its corresponding receiver over a set of parallel 
Gaussian subchannels which may represent time or frequency bins (here we consider transmissions over the 
frequency-selective IC without loss of generality). We denote by Hij{k) the (cross-) channel transfer function 
over the fc-th frequency bin between the secondary transmitter j and the receiver i, while the channel transfer 
function of secondary link i is Hii{k). The transmission strategy of each user (pair) i is the power allocation 
vector Pi = {pi(^)}fcLi over the A^ subcarriers; the power budget of each transmitter i is Ylk=iPii^) — ^i- 
a CR system, additional power constraints limiting the interference radiated by the SUs need to be imposed. 
Here we envisage the use of the following general interference constraints: for each SU i. 



where Wj(/c) G and on £ M!p are nonnegative m-length vectors. Note that constraints in the form of ^ 
are general enough to include, as special cases, for example: i) spectral mask constraints < pi < p™^^, where 
pinax _ (^pf^^^l^k))^^^ is the vector of spectral masks over licensed bands; and ii) interference temperature 
limit-like constraints '^k=i\^pi''^\^)\'^Pii^) — ^pi P ~ where H^'^\k) is the cross-channel 

transfer function over carrier k between the secondary transmitter i and the primary receiver p, and Ipi is the 
maximum level of interference that SU i is allowed to generate. Methods to obtain the interference thresholds 
Ipi when such a knowledge is not available at the SUs' side are discussed in |4H Sec. IV]. For the purposes 



2.1 The SISO case 




(4) 



fc=i 
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of this section, define by 

vr° ^ <! p, G M^v > m{k) <p„ o<Pi< pr^ } , (5) 



= |p. G : Y^Mf^) <P^^ < Pi < pr^l , 



the set of power budget constraint of SU i including exphcitly the power budget and spectral mask constraints. 

Under basic information theoretical assumptions (see, e.g., [HHI)) the maximum achievable rate on link 
i for a specific power allocation profile pi, . . . , p/ is 

where p_.j = (pi, . . . , pj_i, pi+i, . . . ,p/) is the set of all the users power allocation vectors, except the z-th 
one, and + X^jyj \ Hij{^)\^'Pj{^) is the variance of the noise plus the multiuser interference (MUI) over 

subcarrier k measured by the receiver i, with crfik) denoting the power of the thermal noise (possibly including 
the interference generated by the PUs). The system design consists in finding the optimal power allocation 
of the users in order to maximize the information rates of the links, according to some performance metrics. 
The sum-rate maximization [in the absence of the interference constraints (j4])] has been shown to be an NP 
hard problem |42| ; several attempts have been pursued in the literature to deal with the nonconvexity of such 
a problem; however all the proposed schemes are centralized and computationally expensive, which makes 
them non-implementable in a network with no infrastructure, like ad- hoc or CR networks. This motivates 
the formulation of the system design as a NEP: the aim of each player (link) i, given the strategy profile p_, 
of the others, is to choose a feasible power allocation pj that maximizes the rate rj(pj,p_j), i.e., 

maximize rj(pj,p_j) 
Pi 

subject to 

(a): Pi^rr", 1 (7) 

^ I A -p siso 



(b) : Yy^i{k)pi{k) < cXi, 



k=l 



for alH = 1, . . . , /, where 7??^^° and rj(pj, p-i) are defined in ^ and ([6]), respectively. We denote the NEP 
based on (0) by Gsiso = {r^^^°, {ri)i^^) , with r^^^° = Hi^f"^" and Pf^^" being the feasible set of the 
optimization problem ([7]) of SU i. Note that ^siso is an instance of the real player-convex NEP in ([T]). 

Literature review. Special cases of the NEP in ([7]) have been extensively studied in the literature in the 
context of ad-hoc networks, namely when there are only power constraints (a) [H O IH [U |43]. In such a 
simplified setting, given the strategy profile p_j, the optimization problem of each player reduces to: 



maximize '^j(Pi)P-i) 
Pi 



) SlSO 



(8) 



subject to Pi G V{ 

We denote the game resulting from dH) by ^siso = (v^^^°, (ri)f=i), with V^^^° = Ui'^i^^^" ■ Introducing 
hj)lj=i 



the matrices M = (Mi,)f G ]^nixNI ^^^^^ ^ £ defined respectively as 
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\Hu{k)\^ 



2\N 







if i = j; 



Mij = diag 



and [r],,^ 



otherwise, 



(9) 



k=l 



maxfc 



2 ' 



the state-of-the-art-results on ^siso can be collected together in the following theorem, where p{A) denotes 
the spectral radius of A. 

Theorem 1 Given the NEP Qsiso (with no interference constraints), the following hold. 
(ci) Qsiso has a nonempty and compact solution set] 

(b) //M is positive definite, then Qsiso has a unique NE OH]; 

(c) If p{T) < 1, then Qsiso has a unique NE and the Asynchronous Iterative WaterfiUing Algorithm (IWFA) 

based on the waterfilling best-response as proposed in [3] converges to the equilibrium. 

Theorem [1] provides a satisfactory characterization of the NEP ^siso (namely, conditions for the exis- 
tence/uniqueness of the solution and global convergence of distributed solution schemes) under p{T) < 1 
(or M positive definite). However, condition p{T) < 1 may be too restrictive in practice; indeed there are 
channel scenarios resulting in games ^siso having multiple Nash equilibria for which y9(r) <^ 1 does not hold. 
In such cases, the IWFA is no longer guaranteed to converge and there are no algorithms available in the 
literature solving the game ^siso- Moreover, results in Theorem^ as well as the mathematical tools used in 
[H El m El |l3] to study ^siso cannot be applied to the more general t/siso even in the case of unique NE. 
The theoretical analysis of ^siso is then an open problem, which will be addressed in Sec. [71 based on the 
general framework that we introduce in the forthcoming sections. It is worth mentioning that we also derive 
a closed form waterfilling- like expression of the best-response in ([7]) , which is a new results in the literature 
of communications. 

2.2 The MIMO case 

In a MIMO setting, the secondary terminals are equipped with multiple transmit/receive antennas and are 
allowed to transmit over a multidimensional space, whose coordinates may represent time slots, frequency 
bins, and angles. The goal is to find out the most appropriate transmission strategy exploring all available 
degrees of freedom, under alternative interference constraints (in addition to arbitrary power constraints 
imposed on the transmit covariance matrix of the SUs). We envisage the use of interference constraints 
expressed in the following very general form: 

- Null shaping constraints: 



where Qj € C^^^ ^"^j is the transmit covariance matrix of SU i with ny- being the number of transmit 
antennas and Uj G C^^^ ^^^^ is a tall matrix whose columns represent the "directions" along with user 
q is not allowed to transmit. We assume, without loss of generality (w.l.o.g.) that each matrix Uj is 
full-column rank and, to avoid the trivial solution Qj = 0, r{/. < nji- 



UfQi = 0, 
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- Soft and peak power shaping constraints: 

tr (G^Q.Gp,) < i;r and A,nax (F^Q.Fp,) < P = l,2,..., 

which represent a relaxed version of the null constraints by limiting the total average and peak average 
power radiated along the range space of matrix Gpi € C"^'^"'^^ and Fpj G C"^'^"^?', where 1^™ and 
Ip^^^ are the maximum average and average peak power respectively that can be transmitted along the 
directions spanned by Gpi and Fpj. 

Null constraints are enforced to prevent SUs from transmitting over prescribed subspaces, which for 
example can identify: 1) portion of licensed spectrum (the range space of Uj coincides with the subspace 
spanned by a set of Inverse Discrete Fourier Transform (IDFT) vectors); 2) time slots used by the PUs (the 
columns of Uj are a set of canonical vectors); and 3) angular directions identifying the primary receivers as 
observed from the secondary transmitters (the columns of Uj are the set of steering vectors representing the 
directions of the primary receivers as observed from the secondary transmitters). Soft shaping constraints 
can be used instead to control the (average and peak average) power radiated by the SUs along prescribed 
time/frequency/angular "directions" (those spanned by the columns of matrices Gpi and Fpj); for instance, 
classical power constraints, such as per-antenna power constraints [Qjjfcfc < /Sj^ with k = 1, , . . . nx^, or power 
budget constraints tr(Qj) < Pi are example of soft-shaping constraints. In the CR context, null constraints 
are motivated by the interference- avoiding paradigm (also called white-space filling approach): CR nodes 
sense the spatial, temporal, or spectral voids and adjust their transmission strategy to fill in the sensed 
white spaces; such a paradigm has already been adopted as a core platform in emerging wireless access 
standards such as the IEEE 802. 22- Wireless Regional Area Networks (WRANs). Soft-shaping constraints 
find their natural application within the opportunistic communication paradigm, which allows simultaneous 
transmissions between primary and secondary users, provided that the required QoS of the PUs is preserved 
|41| . Under basic information theoretical assumptions (see, e.g., [8]), the maximum information rate on 
secondary link i for a given set of user covariance matrices Qi, . . . , Q/, is 

RiiClu Q-/) = logdet (I + Hf R_i(Q_,)"'HiiQ,) (10) 

where R_,(Q„j) = R„. + ^ HyQjH^ is the covariance matrix of the noise plus MUI, with R„,. G C"^^^"-^* 

denoting the covariance matrix of the thermal Gaussian zero mean noise (possibly including the interference 
generated by the PUs), and assumed to be positive definite; Q_j = (Q^O^yj is the set of all the users covariance 
matrices, except the i-th. one; HjjG C"^!^"^! is the channel matrix between the i-th secondary transmitter 
and the intended receiver, whereas HjjG C"^'^"'^-'' is the cross-channel matrix between secondary source j 
and destination i. Within the above setup, the game theoretical formulation is: for each SU i = 

m^imize i?j(Qj,Q_j) 

subject to 

(a) : tr(Q.O<Pi, 

(b) : Uf Q, = 0, 

(c) : tr (G^Q.Gp,) < K... (f^Q^F^,) < p = 1, 2, . . . , 



fill 



A -^ mirno 



9 



where Qj C !C"'^i^"''^i is an abstract set that can accommodate (possibly) additional constraints on the 
covariance matrix Q^, on top of the power and interference constraints; we only make the (blanket) assumption 
that each Qj is closed and convex. We refer to the NEP based on pT]) as ^mimo = (("P™^™", (Ri)l^^'j, with 
pmimo A Y^^-pmimo ^mimo defined in Note that 0^^^^ is an instance of the complex NEP 

Literature review. The design of MIMO CR systems under different interference-power /interference- 
temperature constraints has been addressed in a number of papers. An overview of centralized solutions 
can be found in |44| , whose approach is to cast the resource allocation of the secondary system into a classical 
convex optimization problem, which is suitable for modeling CR networks composed of one SU only (the case 
of multiple SUs would lead to a NP-hard formulation). Distributed algorithms (mostly) for ad- hoc networks 
based on game theoretical formulations have been proposed instead in \24\ [71 SSI IH [11] ; the state-of-the-art 
result is the asynchronous MIMO IWFA solving the NEP in (jlip . in the presence of constraints (a) [8] and 
(b) |11) only. Results in these papers are strongly based on the specific structure of the optimization problem 
and the resulting solution— the MIMO waterfilling-like expression— and thus are not applicable to the general 
NEP (|lip . ^mimo is thus an other example of a novel game whose solution analysis requires new mathematical 
tools, which is the goal of this paper. The study of ^mimo is addressed in Sec. 17.21 and will result as a direct 
application of the framework developed in the forthcoming sections for complex NEPs. 

2.3 Femtocell OFDM MIMO networks 

Femtocell networks |46) are two-tier systems composed of central macrocell base stations [also termed evolved 
Node-Bs (eNBs) in the 3GPP jargon] and multiple distributed femtocells [also termed as Home evolved 
Node-Bs (HeNBs)]; the latter being short-range, low-power base stations that are deployed and managed by 
the customers at home or in their offices. Due to the massive and dense deployment of HeNBs and their 
uncoordinated nature, interference management is one of the major challenges to be faced in the femtocell 
system design; different interference models and resource allocation algorithms have been proposed in the 
literature (see, e.g., |47| \^8\ [39] and references therein). Here we build on the interference model proposed 
in |49| and introduce a novel spectrum sharing resource allocation and interference management problem in 
OFDM downlink MIMO LTE femtocell systems. Borrowing the idea of CR systems described in the previous 
section, we formulate the design of the radio interference-aware access scheme for HeNBs as a complex player- 
convex NEP, where the eNBs play the role of interference (the PUs) and the HeNBs, acting as SUs, compete 
to optimize their available time/frequency /space degrees of freedom, based on a preliminary sensing and 
prediction of the interference activities of the eNBs. 

Assuming the HeNBs to operate over a frame-base structure (the timeline is divided into T successive 
non-overlapping time-slots), the interference activities of the eNBs over the frequency channels k = 1, . . . ,N 
are modeled as a set of statistically independent first-order homogeneous Discrete Time Markov Chains 
|49| (higher order Markov chain can also be considered without any further conceptual complication); let 
TTf^\t) and 7r^^(t) denote the probability that the subcarrier k in the time-slot t = 1, . . . ,T is idle and busy, 
respectively; introducing the matrix P^*^^ of transition probabilities over channel k, the time evolution of 
the state probabilities Tr^''\t) = [vrf ^(t), 7r^^(t)]^, given the channel state 7r'^'')(0) in the initial time-slot 
t = 0, is 7v^^\t) = 'P^^^7v^''\t — 1). The transition probabilities in the matrix P^'^) can be estimate from the 
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HeNBs either through a knowledge of the traffic statistics of the eNBs, or by sensing the channels and using 
well-known estimation techniques (see, e.g., |50|), or by a combination of these. Hereafter, we assume that 
the HeNBs allocate some of the available slots, let us say the first to time slots of the transmission frame, to 
measure the received interference generated by the active HeNBs and eNBs over each carrier and estimate 
the transition probability matrices P^'^). 

Based on such an estimate and the prediction model 7v^''\t) = P^^'^Tr^''\t — l) of the interference activities, 
the maximum achievable rate on the (downlink) MIMO OFDM channel of HeNB i over T — to consecutive 
time-slots is: 



N 



RiiQi, Q 



T 

E 

t=to + l k=l 



E 



(fc) 
I 



(k) 



(fc) 
B 



(12) 



with Qi 



((Qf''V=i)Lo' 



and r; 



(fc) 



ij 



and r. 



(fc) /p,(fc,t) f^{k,t) 



i,B 



being the achievable rates on 



HeNB link i over carrier k in the time slot t when the channel k is detected as idle and busy, respectively: 



fJJ (Qf'^Q^*)) ^logdet (l + HSf)^R(5,(Qa*))-iH;fQf'*) 



(fc) (^{k,t) f^{k,t)\ A 



'i,B 



Q ,Q 



logdet(l + HPR5^^(Q^ 



v(fc,^)^ 



-ljj(fe)Q(fe,t) 



where R%{Q^^P) ^ E H^J^Qf and RL,_^ 



fik)^ik,t)„ik)H 



(k) 

and R„/g denoting the covariance matrix of the thermal noise and the thermal noise plus the interference 

(k t) 

generated by the eNBs over carrier k at the receiver of the femto-user i, respectively; Q- is the transmit 



(QL'i) = R: 



(fe) 

rii.B 



+ EH;fQf'*)H;f^,withR; 



(fc) 

■riiJ 



J is the cross-channel matrix 



(k) 

covariance matrix of HeNB i over carrier k at time-slot t; and H^^- E 
over carrier k between HeNB source j and the femto receiver i (Hjf^ G C"-f'i''"^i is the channel matrix over 
carrier k of the direct link i). All the channel matrices are assumed to be time-invariant over T consecutive 
time-slots. The goal of each HeNB i is then to maximize the expected MIMO femto-user rate (jl2p . subject to 
power and interference constraints, the latter introduced to preserve the QoS of the users in the macro-cells. 
The proposed game-theoretical formulation is the following: for each HeNB i, 



maximize i?j(Qj,Q- 
subject to 



T N 
t=to+lk=l 



{k,t) 



N 



< Pi. 



k=l 
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N 



^ rave 
— pi 



(k), Xu... (fJ'* )^Qf '^^F^f < i^^r\k), yk, Vt > to, Vp 



k=l 

Qi G Qi, 
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— pi 



N 

E^" 

k=l 



^{k,t)H^{k,t)-^{k,t) 
pi ^i pi 



< vt > to, vp, 



(a) : 

(b) : 

(c) : 

(d) : 

(e) : 

'(13) 

where Qi is an abstract set of possibly additional constraints on Qj = ((Qj-'^'*^)^i)^f(,, assumed to be closed 
and convex. The proposed formulation contains also per-carrier /whole-bandwidth average and peek average 
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power and interference constraints imposed to each HeNB i to protect the users in the macro-cells; the 

(k t) (k t) 

matrices G^j' and F^^' are chosen to select the "directions" along with the power radiated by each HeNB i 
needs to be limited, see Sec. 12.21 for a discussion and interpretation of these constraints. We refer to the NEP 
based on ^ as ^femto = (v^^''^° , {Ri)Zi), with P^^""*" ^ Yl^'Pt""^° and defined in ^femto is 

another example of complex player-convex NEP in ([TJ . 

Literature review. The SISO formulation resulting as special case of (jl3p has been proposed in |49| . in a 
simplified setting where power only budget constraints are considered. Based on [3], the authors in |49) derived 
sufficient conditions for the uniqueness of the NE and the convergence of Jacobi gradient-response algorithms. 
Convergence conditions in |49) are however not easy to be checked (they require an exhaustive search over the 
feasible set of the power constraints) and imply uniqueness of the equilibrium; moreover gradient-response 
schemes suffer from slow convergence (cf . Figure [6]) . The MIMO formulation (jl3p is instead novel and none 
of the approaches proposed in the literature of MIMO games |24| [71 SSI El [11] can be successfully applied 
to design best-response algorithms and study their convergence properties. The analysis of both SISO and 
MIMO formulations can be readily addressed using the framework that we introduce in this paper. 



3 A Theory of (Partitioned) Variational Inequalities 

Our study of NEPs is based on the reduction of a NEP to a VI. The main advantage of this reduction is 
algorithmic, since once the reduction has been carried out we can build on the well-developed VI theory |30| 
in order to design new solution methods for NEPs. A VI representing a NEP obviously has peculiarities 
that must be taken into account; the main purpose of this section then is to present some results on Vis 
geared towards the analysis of NEPs. Special emphasis will be given to Vis defined on a Cartesian product of 
sets, which we term partitioned Vis, as these will be widely used to study NEPs. We note that this section 
contains some results that are new even in the VI literature, such as new sufficient conditions for the block 
P property of a vector function. 

The simplest way to see a VI is as a generalization of the minimum principle for convex optimization 
problems, which is recalled next. Consider a convex optimization problem (in the minimization form), whose 
objective function / : Q i— )■ M is convex and continuously differentiable on the feasible set Q C which is 
a convex and closed subset of M"'. A feasible point x* S Q is an optimal solution of the optimization problem 
if and only if 

(x-x*)^V/(x*) > 0, VxgQ. (14) 

The VI problem is a generalization of the minimum principle (|14p where the gradient V/ is substituted 
by a general real- valued vector function F. More formally, we have the following. Let Q C M'* be a nonempty, 
closed, and convex set and let F : Q — )• M"" be a vector-valued real function. The Variational Inequality VI 
(Q, F) is the problem of finding a feasible point x* G Q such that |30| Def. 1.1.1] 

(x-x*)^F(x*) > 0, VxgQ. (15) 

^When we say that a (vector-valued) function is continuous or continuously differentiable on a closed set we mean that the 
function is so on an open set containing the closed set. 
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The set of solutions to this problem is denoted by SOL(Q, F). As we mentioned earlier, a case that is relevant 
in the analysis of NEPs is that of partitioned Vis. This corresponds to the set Q being a cartesian product 

of lower- dimensional sets: Q = JJ Qj, with each Qi C ]R"» being nonempty, closed, and convex and with 



I 

A 



n = rij. When this structure arises it will be quite natural to partition both F and x accordingly and 
1=1 

therefore write F(x) = (Fj(x))[^;^ and x = (x)^^;^, where Fj : Q — t- M"' is the zth-component block function 
of F and Xj E M"* is the zth-component block of x. 

Several standard problems in nonlinear programming, game theory, and nonlinear analysis can be nat- 
urally formulated as a VI problem; many examples can be found in |30| Ch. 1], [32] . and |33) . Below 
we summarize some known facts and definitions about VI and provide a few new results pertaining to the 
partitioned VI setting. 

3.1 Solution analysis 

In order to present results about the existence and the structure of the solution set of a VI, we introduce 
some function classes. Some of these classes will also play a crucial role in the algorithmic developments in 
later sections. 

Definition 2 A mapping F : Q— t-M", with Q closed and convex, is said to be 

(i) monotone on Q if for all pairs x and y in Q, 

(x-yf (F(x)-F(y)) >0; (16) 

(ii) strictly monotone if for all pairs x ^ y in Q the inequality in ()16p is strict; 

(iii) strongly monotone on Q if there exists a constant Cgm > such that for all pairs x and y in Q, 

(x-yf (F(x)-F(y)) > || x - y f. (17) 
The constant is called strong monotonicity constant. 

If we assume a Cartesian product structure, i.e. F = (Fj(x))^^-^ and Q = Yli Qi, then the function F is said 
to be 

(iv) a Pq function on Q if for all pairs of distinct tuples x = (xj)^^-^ and y = {'yi)l^i in Q, an index i exists 
such that Xj 7^ y^ and 

(x,-yi)^(Fi(x)-F,(y)) > 0; (18) 

(v) a P function on Q if for all pairs of distinct tuples x = (xj)^^-^ and y = (yi)f^]^ in Q, the inequality in 

(fTHI) is strict; 
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(vi) a uniformly P function on Q if there exists a constant CuP > such that for all pairs x = (xj)[^^ and 

y = (yOf=i «^ Q, 

max (xi-yi)^(Fi(x)-F,(y)) > Cup||x-yf. (19) 

1< i <Q 

The constant CuP is called uniformly P constant. 

If a function F enjoys one of the properties above, we will also say that the corresponding VI (Q,F) enjoys 
the property (i.e., ifF is monotone, we say that the VI is monotone, etc.). 

Note that in the case of affine functions, F(x) = Mx + b, there is no difference between strict mono- 
tonicity and strongly monotonicity, and the uniform P property coincides with the P property. Monotonicity 
properties play in the VI realm the same role that convex functions play in optimization. In fact, we recall 
that a differentiable function / is convex, strictly, strongly convex on a convex set Q if and only if its gradi- 
ent is monotone, strictly, strongly monotone on Q. The P properties can be viewed as an extension of the 
monotonicity properties tailored to the possible partitioned structure of the VI; when the partitioned VI has 
only one block, i.e., 1=1, the P properties collapse to the corresponding monotonicity properties. In Figure 
[T]we summarize in a pictorial way some well established relations between these various classes along with 
some of their consequences. Theorem [3] provides instead a formal statement of some existence and uniqueness 
results that will be used throughout the paper. 



strongly monotone = 


> strictly monotone 


=^ monotone = 


> convex solution set 










uniformly P = 


> P 


Po 












unique solution = 


> at most one solution 







Figure 1: Monotonicity and their consequences on Vis. 



Theorem 3 Given the VI (Q, F), suppose that Q is closed and convex and F is continuous on Q. The 
following statements hold: 

(a) The VI(Q, F) has a (possibly empty) closed solution set. If Q is bounded, the solution set is nonempty 

and thus compact f30[ Cor. 2.2.5]; 

(b) //F is monotone on Q, then the VI(Q, F) has a (possibly empty) convex solution set f30\. Th. 2.3.5]; 

(c) IfF is strictly monotone on Q, then the VI(Q,F) has at most one solution f30\. Th. 2.3.3(a)]; the same 

conclusion holds if the VI(Q,F) is partitioned and F is a P function on Q f3(^ Prop. 3.5.10(a)]; 

(d) IfF is strongly monotone on Q, then the VI{Q,F) has a unique solution f30[ Th. 2.3.3(b)]; the same con- 

clusion holds if the VI(Q,F) is partitioned and F is a uniformly P function on Q f30[ Prop. 3.5.10(b)]. 

Note that the uniqueness results stated in part (c)-(d) do not require that the set Q be bounded. For Pq 
partitioned Vis more refined existence/uniqueness results can be obtained; we report them in Appendix lAl 
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3.2 Conditions for the monotonicity and P properties 

The results in the previous subsection and many of the algorithmic developments to follow hinge critically on 
the function F enjoying some monotonicity or P property. However, using directly the definition to establish 
whether such properties hold is in general not possible. It is then useful to be able to establish more practical 
conditions to check the aforementioned properties. It is well known that when Q is an open set and F is 
continuously differentiable on Q, with Jacobian matrix denoted by JF, it holds that |301 Prop. 2.3.2]!^ 



JF(x) >i 0, Vx G Q; 
JF(x) ^ 0, Vx G Q; 
JF-Csml^O, VxG Q; 



(20) 



F(x) is monotone on Q 
F(x) is strictly monotone on Q <;= 
F(x) is strongly monotone on Q 4^ 

where A ^ (A ^ 0) means that A is a positive semidefinite (definite) matrix. The verification of these 
kind of conditions is often difficult and, furthermore, in many practical instances their verification cannot 
easily be linked to physical characteristics of the systems being studied. Therefore, our aim in this subsection 
is that of developing some (conceptually) simpler conditions that permit to deduce the desired F properties 
and that, at least in some instances, can give some further insight into the problem at hand. It is interesting 
to note that, as we will see later on, the conditions we introduce here on pure theoretical grounds also have 
some important algorithmic consequences. 

Let us define the matrix JFiow having the same dimension as JF(x): 



[JF 



low I 



inf [B'^ JF(x)Bl , if r = s, 

— sup I [B^ JF(x) B] ^.^ I , otherwise, 



(21) 



where B € M"^" is an arbitrary nonsingular matrix; and, under the assumption that F and Q have a 
partitioned structure, let us also introduce the "condensed" I x I real matrices Tf and Fp: 



a. 



_ -omax 
-Pij 1 



if i = j, 
otherwise, 



and [Ff] 



with 



inf Aieast(Cf JiFi(x)Ci) and 



max A_ 



0, 



if i = j, 
otherwise, 



sup ||Cf JjFi(x)Cj| 



(22) 



(23) 



where Aieast (A) denotes the smallest eigenvalue of the symmetric part of matrix A, JjFj(x) is the Jacobian 
of Fi(x) with respect to xj, and d € M"'^"' with i = 1, . . . ,/, is a set of arbitrary nonsingular matrices. 
Note that in the definition of Ff we tacitly assumed all a™™ 7^ and P^'^^ are finite; the latter condition is 
equivalent to the boundedness of JjFj(x) on Q. Matrices B and Cj's provide an additional degree of freedom 
in obtaining conditions for monotonicity and P properties of F that can be linked to physical characteristics 
of the systems being studied (see Sec. [7] for some examples). In order to explore the relationship between 
the two matrices Tf and F^t^, we need the following definition. 



■^Conditions in (|20p can be generalized also to the case in which Q is closed; this will be done in Sec. [S] where we introduce 
the VI problem in the complex domain; see Proposition 1291 
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Definition 4 A matrix M G R"^" is called P matrix if every principal minor o/M is positive. 

Many equivalent characterizations for a P matrix can be given. The interested reader is referred to |5H [52] 
for details; here we recall only the following properties, instrumental for our further derivations. Any positive 
definite matrix is a P-matrix, but the reverse does not hold (unless the matrix is symmetric). Furthermore, 
building on the properties of the P-matrices |5H Lemma 13.14], one can show that Tp is a P-matrix if and 
only if p{T-p) < 1, where p(A) denotes the spectral radius of A. 

Matrices JFiow and Tp are useful to obtain sufficient conditions for the monotonicity and P property of 
the mapping F, as given next. 

Proposition 5 Let F : Q— >■ R" be continuously differentiable with bounded derivatives on the closed and 
convex set Q. The following statements hold: 

(a) // JFiow is copositivejfl then F is monotone on Q; 

(b) // JFiow *5 strictly copositive,^ then F is strictly monotone on Q; 

(c) // JFiow is positive definite, then F is strongly monotone on Q with strong monotonicity constant given 

by Csm = Aleast (JFiow) [or Cgm = Aleast (Tf)]0 

If we assume a Cartesian product structure, i.e. F = (Fj(x))f^;^ and Q = Yli Qi, then: 

(d) IfTp is positive semidefinite/PQ-matrix, then F is a monotone/Po function on Q; 

(e) // Tf is a P-matrix [which is equivalent to p{Tf) < 1], then F is a uniformly P-function on Q with 

uniform P constant given by 

CuP(Fj = (24) 

(1 + C(Tf)/5(Tf))2W-^) 

where (^(Tf) — max^^g | [TF]rg|; CL^d 5(Tf) — min{(T([TF]aa) : a C {1,...,/}}, with a{[NL\aa) 
denoting the smallest of the real eigenvalues (if any exists) of the principal submatrix ofM. of order a. 

Remark 6 (On the P-properties) The P property as stated in the proposition above is of a block kind 
(F is partitioned in component blocks of dimension Ui > 1), as opposed to the point kind (F is partitioned 
in component blocks of dimension = 1), as analyzed in [301 Ch. 3]. Furthermore, because of this block P 
property, Proposition [5] does not follow from [301 Proposition 3.5.9] that pertains to the point P property. A 
special case of Proposition [5)^e) can be found in |31j . where the block P property is stated for a mapping F 
associated to some NEPs in the context of CR systems. □ 



matrix A is copositive if x"^Ax > for all x > 0; it is stictly copositive if x-^Ax > for all 7^ x > 0. A positive 
(semi) definite matrix is (strictly) copositive. 

*The least eigenvalue of a real (nonnecessarily symmetric) matrix A, denoted by Aieast(A), is the smallest eigenvalue of the 
symmetric part of A. 
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4 Nash Equilibrium Problems: Fundamentals 



In a standard real NEP there are / players each controlling a variable Xj G M"* that must belong to the 
player's feasible set Qi, which is assumed to be closed and convex: Xj € Qj. In what follows we denote by 
X = (xi,...,x/), the vector of all players' variables, while x_j = (xi, . . . , Xj_i, Xj+i, . . . , x/) denote the 
vector of all players' strategies variables except that of player i. The aim of player i, given the other players' 
strategies x_j, is to choose an Xj G Qi that minimizes his cost function /j(xi,x_j), i.e., 

minimize /i(xj, x_j) 

(25) 

subject to Xj G Qi. 

Note that the players' optimization problem are coupled since the players' objective function (may) depend on 
the other players' choices. Define the joint strategy set of the NEP by Q = ^[=1 2i> whereas = Yijj^i Q-ji 
and set f — (/«)[=!■ The NEP is formally defined by the tuple Q = {Q, f). A solution of the NEP is the well- 
known Nash Equilibrium (NE), which is formally defined in ([2]). 

A useful way to see a NE is as a fixed-point of the best-response mapping for each player; this suggests 
the use of (iterative) best-response-based algorithms to solve the game. Given the limitations of classical 
fixed-point results in the study of convergence of best-response based algorithms (cf. Sec. [1]), we address this 
issue by reducing a NEP to a VI problem, and then using the well developed theory of Vis |30| . This approach 
has the advantage of leading quite naturally to the derivation of implementable distributed algorithms along 
with their convergence properties. 

4.1 Connection to variational inequalities 

At the basis of the VI approach to NEPs there is an easy equivalence between a real NEP and a suitably 
defined partitioned VI. This equivalence follows readily from the minimum principle for convex problems and 
the Cartesian structure of the joint strategy set Q |30| Prop. 1.4.2]. 

Proposition 7 Given the real NEP Q = (Q, f), suppose that for each player i the following hold: 

i) the (nonempty) strategy set Qi is closed and convex; 

ii) the payoff function /j(xj,x_j) is convex and continuously differentiable in Xj for every fixed x_j. 
Then, the game Q is equivalent to the VI(Q, F), where F(x) = (Vxj/i(x))[^-i^. 

In the sequel we refer to the VI(Q, F) defined in previous proposition as the VI associated to the NEP Q. 
It is possible to relax the assumptions in Proposition [7] and still get useful connections between games and 
Vis |20| : but since our aims are mainly computational, we don't pursue this topic further. Indeed, throughout 
the paper, we will make the following blanket convexity/smoothness assumptions, unless stated otherwise. 

Assumption 1. For each i = 1, . . . ,/, the set Qi is a nonempty, closed, and convex subset of M"' and the 
function fi{xi,x^i) is continuously differentiable on Q = Y\- Qi and convex in Xj for every fixed x_i G Q~i- 

Assumption 2. For each i = 1, . . . , I, each function /i(x) is twice continuously differentiable with bounded 
derivatives on Q = J^- Qi. 
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4.2 Existence and uniqueness of a NE 

Building on the VI reformulation in the previous section and the existence/uniqueness results for partitioned 
Vis given in Theorem [3l we can easily state the following theorem that needs no further proof. 

Theorem 8 Given the real NEP Q = (Q,f), suppose that Q satisfies Assumption 1 andletF{x) = (Vxj/i(x))[_ 
Then, the following statements hold: 

(a) Suppose that for every i the strategy set Qi is bounded. Then the NEP has a nonempty and compact 

solution set; 

(b) Suppose that F(x) is a monotone function on Q. Then the NEP has a convex (possibly empty) solution 

set; 

(c) Suppose that F(x) is a P (or strictly monotone) function on Q. Then the NEP has at most one solution; 

(d) Suppose that F(x) is a uniformly-P (or strongly monotone) function on Q. Then the NEP has a unique 

solution. 

Sufficient conditions for F(x) being a strictly, strongly monotone or a (uniformly- )P function on Q are 
given in Proposition [5j The condition given in Theorem Ufa) guaranteeing the existence of a solution for a 
NEP requires that the strategy set of each player be bounded. More general results on the existence of a 
solution for a NEP with unbounded strategy sets can be found in Appendix [Al Note that the (existence and) 
uniqueness results stated in parts (b)-(d) do not require that the set Q be bounded. 

Remark 9 (On the uniqueness conditions) According to Proposition [5{b), and supposing that F(x) = 
(Vx,/j(x))[^^ is continuously differentiable with bounded derivatives on Q (Assumption 2), a sufficient con- 
dition for the uniqueness of the NE is that the matrix Tf defined in (j22p is a P matrix. It turns out that 
this condition is sufficient also for global convergence of best-response asynchronous distributed algorithms 
described in Sec. [5l Note that if Tf is a P matrix, it must be a™'° = infzgQ ['^mm(Vx-/j(z))] > for all i, 
where Aniin(Vx^/j(z)) denotes the minimum eigenvalue of Vx-/i(z). Thus an implicit consequence of the P 
assumption on the matrix Tf is the uniform positive definiteness of the matrices V^j/i on Q, which implies 
the uniformly strong convexity of /j(-, x_j) for any given x_i E Q_i and thus the uniqueness of the solution 
of the i-th player's optimization problem, for any given x_i G Q-i. This latter implication is a quite natural 
consequence of the uniqueness of the solution of the NE: one cannot expect the game to have only one NE if 
the players' problems have multiple solutions. 

The observations above also pave the way for a nice interpretation of the P condition. Note first that the 
;5's in the definition of the matrix Tf measure the coupling of the players' optimization problems: the larger 
the /3's, the more coupled the players' subproblems are. Indeed, if all the /3's were 0, the game Q =<Q, f> 
would decompose into / uncoupled optimization problems; in such a case, requiring the matrix Tf to be P 
simply amounts to requiring all a's to be positive, which obviously implies uniqueness of the solution. It is 
reasonable that if the /S's increase from zero but remain small enough with respect to the a's, the game will 
still have a unique solution. The P property quantifies how large the /?'s can grow while still preserving the 
uniqueness of the solution. □ 

We conclude this subsection providing a sufficient condition for the matrix Tf in (f22|) to be a P (positive 
definite) matrix, which can be derived by elementary diagonal dominance arguments. 
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Proposition 10 The matrix in (|22p is a P-matrix if one of the following two sets of conditions are 
satisfied: for some w = {wi)l^^ > 0, 

Wi ^ a™" ui,- ^-^ a™™ 

If actually both conditions in (j26p are satisfied, then Tf «s positive definite. 

The sufficient conditions in Proposition [10] will be sfiown in Sec. [7] to have an interesting physical interpre- 
tation in the context of power control problems in CR systems. 

4.3 Problem classes 

Based on the previous results, it is natural to introduce the following classes of real NEPs, which will play a 
fundamental role in the rest of the paper. 

Definition 11 A real NEP G = (Q,f) is: 

i) a monotone Nash equilibrium problem if Assumption 1 holds and the associated 1//(Q, F) is monotone; 

ii) a uniformly P Nash equilibrium problem if Assumption 1 holds and the associated VI{Q,Y) is uniformly 

P; 

iii) a Pf Nash equilibrium problem if Assumptions 1 and 2 hold and the matrix Tf of the associated 

VI{Q,¥) IS P. 

From the development so far, Figure [2] summarizing the relations between these classes of problems needs no 
further explanation, except for the fact that P^ (and thus uniformly P) NEPs are not a subclass of monotone 
NEPs. This is shown by the following example. 




Figure 2: Relation among NEP classes. 



Example 12 (A P-f NEP which is not monotone) Consider a real NEP with two players, each con- 
trolling one scalar variable: xi and X2. The players' problems are 

■ ■ ■ 1 2 . ■ • • 1 2 1 

mmimize^ji — + 4x1X2 mmimize^jj —X2 — —3^1X2 

subjectto xi G [0, 10] subjectto X2 G [—2, 2] 
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The VI associated to this is VI([0, 10] x [-2,2],F), with F = [xi +Ax2,X2 - (l/8)xi]'^. The symmetric part 
of JF, JF^, and the matrix Tf are given by: 

1 -4 " 

-1/8 1 

Since JF^ has a negative determinant, F cannot be monotone; on the other hand it is easy to check that the 
two principal minors of Yp are positive, implying that ^ is a P NEP. 

By Theorem [8l a uniformly P (and thus a Pt) NEP has a unique solution, whereas monotone NEPs have 
a convex solution set that can be empty. While obviously a NEP with a unique solution is not necessarily a 
P NEP, the latter class of problems covers an interesting array of problems, as we will be shown later on. 

Centralized algorithms for monotone and uniformly P-f NEPs, based on VI theory, are well-known vol 
II]; in this paper, we focus on the more challenging issue of devising distributed (and possibly asynchronous) 
solution schemes for NEPs, which is the topic of the next section. 

5 Distributed Algorithms for NEPs 

This section along with the next one constitute the core theoretical part of the paper. We develop here a 
novel theory that allows devising distributed algorithms for computing Nash equilibria in several significant 
settings. More specifically, we will provide novel distributed (asynchronous) algorithms for the solution of 

(a) monotone NEPs; 

(b) Pt NEPs. 

Since monotone NEPs may have multiple solutions, in case (a), we will further consider both the situations 
in which one is interested in computing any one solution, and the situations in which one wants to select 
the best solution, according to a given criterion. Settings (a) and (b) cover a wide range of situations and 
they will be shown in Sec. [7] to encompass several important equilibrium models in signal processing and 
communications; in particular the models proposed in Sec. [2j In each of the settings above we will provide 
best-response-based distributed algorithms along with their convergence properties; the proposed algorithms 
differ in: i) the computational effort; ii) the players' synchronization/signaling requirements; and iii) the 
convergence speed. Note that while centralized solution methods are known for uniformly P NEPs, the 
development of distributed algorithms for this class of games is at the time of this writing an open problem. 

This section is organized in three parts. Sections 15.11 and 15.21 focus on algorithms for and monotone 
NEPs, respectively; results in this sections will be the building blocks for the more difficult issue of equilibrium 
selection problem addressed in Sec. 15.31 

5.1 Best-response distributed algorithms for P-f NEPs 

Since in a NEP every player is trying to minimize his own objective function, a natural approach to compute 
a solution of a NEP is to consider an iterative algorithm wherein all the players, given the strategies of the 
others and according to a given scheduling (e.g., sequentially or simultaneously), update their own strategy 



JF. 



1 31/16 
31/16 1 
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by solving their optimization problem ()25p . Here, we focus on a very general class of best-response-based 
algorithms, namely the totally asynchronous best-response algorithms (in the sense specified in [53j). In these 
schemes, some players may update their strategies more frequently than others and they may even use an 
outdated information about the strategy profile used by the others; which is very appealing in many practical 
multiuser communication systems, such as wireless ad- hoc networks or CR systems wherein synchronization 
requirements are hard to enforce. 

To provide a formal description of the algorithm, we need to introduce some preliminary definitions. In 
an asynchronous scheme, the users may not update their own strategies at each iteration; let denote then 
by 7i ^ T ^ {0, 1,2, . . .} the set of times at which player i updates his own strategy Xj, denoted by x^"^ 
(thus, implying that, at time 7i ^ 71, x^"^ is left unchanged). Moreover, in computing their optimal strategy, 
the users can use an outdated version of the others' strategies; let then r* (n) be the most recent time at 
which the strategy profile of player j is perceived by player i at the n-th iteration (observe that T^j{n) satisfies 
< T* (n) < n). Hence, if player i updates its strategy at the n-th iteration, then he minimizes his cost 
function using the following (possibly) outdated strategy profile of the other players: 

A^Kn)) A ({r\{n)) (rU,(n)) (r^+iW) (r]{n))\ 

Some standard conditions in asynchronous convergence theory, which are fulfilled in any practical implemen- 
tation, need to be satisfied by the schedule Ti's and r* (7i)'s, namely for each i: 

Al) < T^j{n) < n (at any given iteration n, each player i can use only the strategy profile x^. ^"^^ adopted 
by the other players in the previous iterations ); 

A2) limfc^oo ''"j(^fc) = +00, where {n^,} is a sequence of elements in % that tends to infinity (for any given 
iteration index n^, the values of the components of x_- in (j27p generated prior to are not used 



in the updates of x^"'^ when n becomes sufficiently larger than n^); 
A3) |7i| = oo (no player fails to update his own strategy as time n goes on). 

Using the above definitions, the totally asynchronous algorithm based on the best-responses of the players is 
described in Algorithm [TJ The convergence properties of the algorithm are given in Theorem [T3j 

Algorithm 1: Asynchronous Best-Response Algorithm 

(5.0) : Choose any feasible x^*^) G Q and set n = 0. 

(5.1) : If x^") satisfies a suitable termination criterion: STOP 
(S . 2) : for ? = 1, . . . , /, compute 



x^ G argmin/i (xj, x^. ^"^M , lin^Ti 
x,-"^ , otherwise 



(28) 



(S.3) : n ^ n-M; go to (S.l). 
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Theorem 13 Let Q = {Q, f) be a P-f NEP. Any sequence {x^^^jJ^Q generated by AlgorithmUl converges to 
the unique NE of Q, for any given updating schedule of the players satisfying assumptions A1-A3. 

Proof. See Appendix [Bj ■ 

Remark 14 (Flexibility of the algorithm) Algorithm [1] contains as special cases a large number of algo- 
rithms, each one obtained by a possible choice of the schedule of the users in the updating procedure (i.e., the 
parameters {T*(n)} and {Ti}). Examples are the simultaneous (Jacobi scheme) and sequential (Gauss-Seidel 
scheme) updates, where the players update their own strategies simultaneously and sequentially, respectively. 
Indeed, the Jacobi update corresponds to the schedule T^j{n) = n and Ti = {1,2, . . .} for all i and j, whereas 
the Gauss-Seidel scheme is obtained by taking T* (n) = n and T = {i, i + 1, i + 21, . . .} for all i and j. More- 
over, variations of such a totally asynchronous scheme, e.g., including constraints on the maximum tolerable 
delay in the updating and on the use of the outdated information (which leads to the so-called partially 
asynchronous algorithms), can also be considered |53| . An important result stated in Theorem [13] is that all 
the algorithms resulting as special cases of Algorithm [1] are guaranteed to reach the unique NE of the NEP, 
under the same set of convergence conditions, since the matrix Tf does not depend on the particular choice 
of {Tj(n)} and {T}- Note that all the algorithms coming from Algorithm [1] are robust against missing or 
outdated updates of the players. This feature strongly relaxes the constraints on the synchronization of the 
players' updates; which makes this class of algorithms appealing in many practical distributed systems. 

Note that the (synchronous) projection-response algorithms for monotone Vis (and thus NEPs) proposed 
in [31] and [301 133] are not guaranteed to converge if applied to a P-r NEP that is not monotone. □ 

Remark 15 (On the convergence conditions of best-response algorithms) Global convergence of Al- 
gorithm [1] is guaranteed under the P property of Tf (or equivalently p{Ty) < 1). We have already pointed 
out however that such a condition cannot be satisfied if there is a player whose cost function has a singular 
Hessian, even in just one point. In fact, if this is the case, we have, a™™ = for some i, let us say i = 1 
without loss of generality, which implies that the matrix Ff has a 1 in the left-upper corner. Since the matrix 
Fp is nonnegative, we have that this implies /o(Ff) > 1 Th. 1.7.. 4]. Assuming that the element 1 is 
contained in an irreducible principal matrix, we will actually have p{T-p) > 1. Note that the irreducibil- 
ity assumption is extremely weak and trivially satisfied if the matrix Ff is positive, which is true in many 
applications. In the next section we discuss a remedy for this issue. □ 

5.2 Proximal distributed algorithms for monotone NEPs 

In this section we deal with monotone NEPs (see Definition [TT]) . The main difference between Px (or 
uniformly P) NEPs and monotone NPEs is that the former always have a unique solution, while monotone 
NEPs may have multiple solutions, and the players' objective functions need not be uniformly strongly 
convex; see Remark 1151 It turns out that, for such games. Algorithm [1] may fail to converge. There is a host 
of solution methods available in the literature to solve monotone real Vis and thus monotone NEPs (see, e.g., 
|3m Vol. II]), but these algorithms are centralized. Recently, in |34| . the authors proposed some distributed 
synchronous schemes for solving monotone Vis, based on the gradient-response mapping; we have already 
discussed the main drawbacks of these algorithms, see Sec. [1] (see also Sec. [7] for some numerical results). 
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The development of distributed best-response algorithms for solving monotone NEPs with (possibly) 
multiple solutions is a challenging task; in this subsection, we provide a first answer to this issue building 
on a regularization technique known as proximal algorithms, see [301 Ch 12] for an introduction to proximal 
point methods for Vis. The proposed approach is to reduce the solution of a single monotone NEP to the 
solution of a sequence of Pf NEPs with a particular structure. The advantage of this method is that we 
can efficiently solve each of the P-f NEPs with convergence guarantee using Algorithm 1 introduced in the 
previous section; the disadvantage is that, to recover the solution of the original single monotone NEP, we 
have to solve a (possibly infinite) number of P-f NEPs. However, it is important to remark from the outset 
that this potential drawback is greatly mitigated by the fact that, as we discuss shortly, (i) we only need 
to solve these Pf NEPs inaccurately; (ii) the (inaccurate) solution of the Px NEPs usually requires little 
computational effort; and (iii) in practice, a fairly accurate solution of the original NEP is obtained after the 
solution of a limited number of Pf NEPs. 

Before introducing the formal description of the algorithm, let us begin with some simple observations 
motivating how the sequence of Pf NEPs is built. Let ^ = (Q, f ) be a monotone NEP; consider a perturbation 
of this game defined as Q-r^y = (Q, {fi + (t/2)|| • — yi|P)f=i), where r is a positive parameter and y = (yj)f=i 
is a given vector in M", with each yj E M"*; we term y center of the regularization. Note that QT-,y is the 
game wherein each player i, anticipating x_j E Q-i, solves the following convex optimization problem: 

minimize /j (xj , x_ j ) + ^ 1 1 Xj - y^ | p 

(29) 

subject to Xj E Qj. 

Let us consider now the VI reformulations of Q and GT,y, given by VI(Q, F) and VI(Q,FT-^y) respectively, 
where F(x) = (Vx;/j(x))^^-^ and FT-^y(x) = F(x)+r (x — y), and let us introduce the matrices Tp and Tp^^ 
associated to Q and GT,y, respectively. It is not difficult to check that 



■Ft,- 



Tp + rl. (30) 



Note that Tp^ ^ does not depend on y. It follows readily from (j30p that if r is large enough, Tp^ ^ is a P 
matrix, meaning that Gr,y is a Pf NEP, for any given y E M". More specifically, using the definitions of 
/S^'^^'s and a™™'s as given in ([23|) . we have the following. 

Lemma 16 Let Q = (Q, f) be a monotone NEP. For any given y E M", the game Q^-^y = (Q, (/j + (t/2)- 
II • ~yj|P)f=i) P-f NEP for every r larger than f (independent ofy), with 



T = max 

l<i<I 



5^/3---ar^. (31) 



Although the assumption that Q be monotone will play a key role in the forthcoming derivations, it is 
important to remark that what we need for the above lemma to hold is only that the cost functions /j's in 
the game have bounded second order derivatives on Q. 

Nice as it is, the result above would be of no practical interest if we were not able to connect the solutions 
of ^T- y to those of Q. Indeed, the solutions of Q and QT-,y are in general different but, nevertheless, there 
exists a connection between them: a point x* is a solution of Q if and only if x* is a solution of Gr,:x.*- 
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Proposition 17 Let Q = {Q, f) be a monotone NEP. For any given t > 0, £ Q is a solution of Q if and 
only if X* is a solution ofQr^yi*. 

Proof. See Appendix O ■ 

Lemma [T6] and Proposition [T7] open the way to the design of convergent distributed algorithms for monotone 
NEPs, as shown next. Let us choose r being large enough so that Qr^y is a P-f NEP (cf. Lemma [T6]) . It 
follows from Theorem [8] that Qr^y has a unique solution, denoted by ST-(y)- Using S,-(y), Proposition 1171 can 
be restated as follows: x* is a solution of Q if and only if it is a fixed point of Sr(»), i.e., x* = Sr(x*). It seems 
then natural to compute the solutions of Q using the fixed-point-type iteration x^""*"^) = St-(x(")), starting 
from a feasible point x^"-*; which corresponds to solving the sequence of NEPs Q^y.{n) for n = 0, 1, . . .. If r is 
sufficiently large [e.g., as in (I3ip ]. each Q^^{n) is a P-f NEP (cf. Lemma [T6]l . and thus its unique solution can 
be computed in a distributed way with convergence guarantee by the asynchronous best-response algorithm 
described in Algorithm [1] (cf . Theorem [T3|) . The above discussion motivates the following algorithm for 
computing the solutions of a monotone NEP, whose convergence properties are given in Theorem [18] below. 

Algorithm 2: Proximal Decomposition Algorithm (PDA) 
Data : Let r > be given. 

(5.0) : Choose any feasible x^''^ G Q and set n = 0. 

(5.1) : If x^"^ satisfies a suitable termination criterion: STOP. 

(5.2) : Solve the game g^,,{„) and set x("+i) = S^(x(")) 

(5.3) : n ^ n -M; go to (S.l). 

Theorem 18 Let Q = {Q, f) be a monotone NEP with a nonempty solution set. Suppose that r is large 
enough so that Tf^ ^ is a P matrix. Then, Algorithm\^is well defined, and the sequence {x^^^jJ^Q generated 
by the algorithm converges to a solution of the game Q. 

Algorithm [2] is of great conceptual interest, but its applicability is limited, unless one is able to easily 
compute St-(x(")). Although there are interesting problems in which this can be done efficiently (see Sec. 
[7|), in general one is expected to solve a number of P-j- NEPs, each of them requiring an infinite iterative 
method to compute each St-(x(")). To overcome this issue, we propose next a variant of Algorithm [21 in 
which suitable approximations of Sr(x^")) can be used. Algorithm [3] below describes such a variant, where 
we have added a further degree of freedom in the updating rule: the new iteration x^""*"^) is not necessarily 
given by (an approximation of) St-(x^")), but lies instead on the line connecting the old iteration x*^*^) to (the 
approximation of) St-(x''")). 

Algorithm 3: Approximate Proximal Decomposition Algorithm (APDA) 

Data : Let {e^^^j^^g, {ry^^^^^o and r > be given. 

(5.0) : Choose any feasible x^*'^ G Q and set n = 0. 

(5.1) : If x^"^ satisfies a suitable termination criterion: STOP. 

(5.2) : Solve the game Q^^^(r.) within the accuracy e("): Find a z(") s.t. ||z(") - S^(x('"))|| < e("). 

(5.3) : Set x("+i) = (1 - r/("))x(") +r?WzW. 

(5.4) : n ^ n + 1; go to (S.l). 
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The error term e*^") measures the accuracy used at iteration n in computing the solution St-(x("')) of G^- x(") • 
The parameter j/"-* instead, introduces a memory in the updating rule, it establishes where exactly we move 
along the line passing through the old iterations x^") and z^") . Note that if we take e^"^ = and 7/^"^ = 1 for 
all n, Algorithm [3] reduces to Algorithmic The advantage of Algorithm [3] with respect to Algorithm [2] is that 
z^") can be computed in a finite number of steps, so that Algorithm |3] becomes implementable in practice. 
Obviously, the errors e(")'s and the parameters r/^'^^'s must be chosen according to some suitable conditions, 
if one wants to guarantee convergence. These conditions are established in the following theorem. 

Theorem 19 Let Q = {Q, f) be a monotone NEP with a nonempty solution set. Suppose that r is large 
enough so that Tf^ y ^■^ ^ P-matrix. Choose {e^"-*} C [0, oo) such that Yln!=i ^''"^ ^ ^ '^'^^ {r/*-"^} C [Rrm Rm], 
with < Rm < Rm < 2. Then, Algorithmic is well defined, and the sequence {x^"^}5^q generated by the 
algorithm converges to a solution of Q. 

The proof of Theorem [19] (and thus also Theorem [TH]) is a consequence of the following facts and thus 
is omitted: i) |30| Th. 12. .3. 9]; ii) The observation that Q is equivalent to the VI(Q, F) (Proposition [7| , 
with F = (Vx,/j)f=i, and the VI(Q, F) has a solution; and iii) Under the P property of Tf^ y. Step S.2 of 
Algorithm [3] (and Algorithm [2]) is well defined (Theorem [13]). 

It is interesting to remark that for sake of simplicity we assumed r to be a fixed number. However, r can 
be varied from iteration to iteration provided that r G (r, r™*^^], where r™'^^ is any finite number. 

While the utility of having the possibility to use inexact solutions in Step S.2 of Algorithm [3] is apparent, 
the usefulness of Step S.3 is less evident. This kind of "averaging" is known as over-relaxation and has its 
roots in classical successive over-relaxation methods for solving systems of linear equations [551 Sec. 7.4]. In 
our context, the extra degree of freedom offered by Step S.3 can bring numerical improvements; see, e.g., |30) . 

Algorithm [3] is conceptually a double loop scheme wherein at each (outer) iteration n, given x^"'^ one 
needs to compute the approximation z^"^ which requires an inner iterative process. Since the condition 
e^"^ < oo implies e^'^'^ ^ 0, when the iterations progress, Sr(x^'^'') has to be estimated with an increasing 
accuracy. However, in practice, this is not a problem since when iterations progress, {x*-"^} usually converges, 
implying Hx^""''-'^) — x(")|| — )• 0. One can then use x^^^^ as a good approximation to initialize any (inner) 
procedure in Step 2 to compute z^"^. It turns out that, in spite of the increasing precision requirements, in 
practice a suitable z^"^ in Step 2 can be computed very easily. 

Finally, observe that a natural choice for computing z*^"^ in Step 2 of Algorithm [3] is Algorithm [1] When 
this choice is made. Algorithm [3] also becomes an asynchronous method, having all the desired features 
described in the previous section. The only difference with Algorithm [1] is that, in Algorithm [3] "from time 
to time" (more precisely when the inner termination test ||z^"^ — St-(x^"'')|| < e^"-* in Step 2 is satisfied) 
the objective function of the players are changed by updating the regularizing term from ^||xj — x-"^p to 
^||xj — x^""^^^ P, which generally requires some coordination among the players to establish when a satisfactory 
approximation z^") has been reached. The remark below deals exactly with issues related to this aspect. 

Remark 20 (On the Inner Termination Criterion) We have seen that, in Step 2 of Algorithm [3] the 
players must be able to decide whether ||z(") — St-(x("))|| < e^"^ holds. In the following we suggest a simple 
distributed protocol to do that. Observe preliminarily that an error bound on the distance of the current 
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iteration z^") from the solution St-(x*^"')) of ^f") '^^^ obtained by solving a convex (quadratic) problem 
(see, e.g., |30| Prop. 6.3.1], [30| Prop. 6.3.7]). For example, under the P properties of Tp defined in 
(j30p . the following error bound holds for the game G^^(n) \30\ Prop. 6.3.1]: a (finite) constant c > exists 
such that 

||z-S,(x("))|| <c||F;"*(z)|| VzeQ, (32) 

where F"^^^ (z) = z — Ilg (z — F(z) — rz), with IIq (x) denoting the Euclidean projection of x onto the 
closed and convex set Q. Note that, since Q has a Cartesian structure, F"^* (z) can be partitioned as 
pnat 1^2) = (I^F"*^* (z)] .)^_-^, where each ^F^^^ (z)]^ = Zj — IIq. (zj — Fj(z) — rzj) can be locally computed by 
the associated player i by solving a quadratic programming, as long as Fj(z) is available at the player side. 

Using ()32p . the implementation of the Step 2 of Algorithm [3] can be obtained as follows. Each player i 
choses preliminarily a suitable local termination sequence {e^^^jn C [0, 00) such that X^^x e^"^ < 00; the 
termination criterion of each player i becomes then || [F"^* (x("))]J| < e-"\ which can be locally imple- 
mented. Once the desired local accuracy is reached by all the players, they can all update the center of their 
regularization. Note that this protocol guarantees that the requirement on the sequence e^*^-* in Step 2 as 
stated in Theorem 1191 is met, since e^") = satisfies X^^i e^*^^ < 00. The last issue to address for 

a practical implementation of this protocol is to understand how the players can detect the others having 
reached the desired termination criterion. This can be done by exchanging one bit of information, if some 
signaling is allowed; otherwise each user can just update its regularization after experiencing no changes in 
11 i^pnat (x^"))] Jl (or his cost function) for a prescribed number of iterations. □ 

5.3 Equilibrium selection for monotone NEPs 

In the previous section we discussed distributed algorithms for the computation of a solution of monotone 
NEPs. A feature of these algorithms is that they converge under mild conditions that do not imply the 
uniqueness of the NE of the NEP. In the presence of multiple equilibria, however, the proposed algorithms do 
not allow to perform any selection of the solution they reach, but they may converge in principle to any NE 
of the game; which makes the achievable system performance unpredictable. It would be interesting instead 
to be able to select, among all the solutions of the game, the one(s) that satisfies some additional criterion. 
We refer to this problem as equilibrium selection problem. In this section, we address this issue; the outcome 
will be a novel set of distributed algorithms along with their convergence properties that solve the equilibrium 
selection problem; this additional feature comes at the price of a (moderate) increase of the complexity in 
computing the players' best-response solution and signaling among the players. 

Let us introduce first an informal description of the algorithm. Let Q =< Q, f > be a monotone NEP and 
let SOL(Q, f) denote its solution set, assumed to be nonempty without loss of generality. Recall that since 
Q is monotone, SOL(Q, f) is always convex (cf. Theorem [8]). Stated in mathematical terms, the equilibrium 
selection problem consists in solving the following bi-level optimization problem: 

^In many practical applications, as those considered in the second part of the paper, each [F"''' (z)]^ can be computed by 
local measurements from the players. For instance, this happens when the strategies of the other players affect the cost function 
of each player by a term that is locally measurable as additive interference noise. 
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minimize <Pi^) 

X 

subject to X G SOL(Q, f), 



(33) 



where the function (p : M" — )■ M is assumed to be continuously differentiable and convex. The function (p 
thus defines the additional criterion according to which one wants to select a solution in the set of the Nash 
equilibria of the game Q: solving (j33p indeed corresponds to choosing the NE of Q that minimizes (p. Note 
that, under the monotonicity of Q, (j33p is a convex optimization problem [SOL(Q, f) is convex]. However, 
standard solution techniques cannot be applied because the feasible set SOL(/C, f ) is only implicitly defined 
and, in general, it is not expressed as a standard system of inequalities. To overcome this difficulty, and in the 
same spirit of the previous section, instead of attacking problem (|33p directly, we propose to solve a sequence 
of standard regularized NEPs ("standard" means a game whose players' feasible sets are not of an implicit 
type, and therefore can be solved by classic methods, like Algorithm [1]) . Each standard regularized game 
has the following structure GT,e,y = (Q, (/« +£(/)+ {t/2)\\ • — yi|P)f=i), where e and r are fixed positive 
constants and y = (yj)f=x ^ given point in M" with each E R"'; ^T,e,y is a NEP wherein each player i, 
anticipating x_j G Q-i, solves the following convex optimization problem: 

minimize fi{xi, x_j) + e(^(xj, x_i) + ^Hxj - yj|p 

(34) 

subject to Xj € Qi. 

Note that the players' problems in this game differ from (j29p in the presence of the additional term £(j){xi, x_j) 
in the objective function. It is not surprising that the function (j) appears in the objective function of the 
players; roughly speaking, it represents the additional amount of information to be included in the game to 
"drive" the system toward the desired solution. 

Proceeding as in the previous section, we can now establish the connection between the regularized NEPs 
QT,e,y 3.nd the equilibrium selection problem (133p . Fist of all note that, in the setting of problem (|33l) . the 
NEPs ^T,£,y is equivalent to the VI(Q,FT- £ y) where F^^^ y = F + e V(/) + t (I — y) and I : x i— t- x is the identity 
map. Denoting by T^^^y the matrix defined in (f22|) and associated to F^-^^^y, Lemma [21] below shows that 
there exists a sufficiently large r such that Tp^ ^ ^ is a P matrix, implying that ^r,e,y is a Pf NEP. 

Lemma 21 Let Q = {Q, f) be a monotone NEP; let </> : i— )• M 6e a continuously differentiable function on 
Q whose gradient V(/) is Lipschitz continuous on Q, with constant L^; and let e > be given. For any fixed 
y G M"-, the game ^r,£,y = (Q, (fi + (t/2)|| • -yi|P)f=i) with e £ [0,e] is a Pr NEP for every r larger 
than Tg {independent on y and e) 

f, ^ max jZ/^r^ - «r I + U - (35) 

where fSf^'''' 's and af"" 's as defined m (Eg). 

Under the setting of Lemma [2T| the game ^T,e,y is a Px NEP and thus has a unique solution, denoted 
by ^r,e (y ) [cf • Theorem [8] ; such a unique solution ^r,e (y ) can be computed with convergence guaranteed 
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using Algorithm [T] on GT,e,y The solution of the original equilibrium selection problem (j33p can be recovered 
using the fixed-point-type iteration x*-"^"*^^ = ^(n) (x^"-*), starting from a feasible point x*-^-* and by suitably 
varying e^"); which corresponds to solve the sequence of NEPs ^(„) ^(„) for n = 0, 1, . . .. This procedure is 
made formal in Algorithm 3] below. 

Algorithm 4: Proximal-Tikhonov Regularization Algorithm (PTRA) 

Data : Let {e^")} | and r > be given. 

(S.O) : Choose any feasible x^'^^ G Q and set n = 0. 

(S. 1) : If x*-"-* satisfies a suitable termination criterion, STOP. 

(5.2) : Set x^'"+^^ to be the solution of the game ^(n) ^l")- 

(5.3) : Set n+1 and return to (S. 1). 

Note that, since in Algorithm [4] the sequence {e^"^} converges to zero, there always exists an e > such 
that e^"^ G [0,e], implying by Lemma [2T] that for a sufficiently large r all the games e(") x(") Step 2 

of the algorithm are Py NEPs and thus have a unique solution; this makes the sequence {x(")}J^g generated 
by Algorithm |4] well defined. The convergence properties of the algorithm are given in the following theorem. 

Theorem 22 Let G = {Q, f) be a monotone NEP with a nonempty solution set SOL(Q, f). Consider the 
equilibrium selection problem (j33p and suppose that i) (j) is continuously differentiate and convex on Q; ii) 
the level sets of (j) on SOL(Q, f) are bounded; and Hi) V(/> is Lipschitz continuous on Q, with constant L^. 
Moreover, suppose that t is large enough so that is a Pf NEP for any n, and choose the sequence 

{e^*^)} such that e^") > for all n, \, 0, and Yl^=o ^^^'^ — Then Algorithm^ is well defined; the 

sequence {x^^^jJ^Q is bounded; and every of its limit points is a solution of (j33p . 

Proof. See Appendix [Dj ■ 

Theorem [22] guarantees convergence of the algorithm under mild assumptions. Conditions on (p are pretty 
standard; in particular, assumptions i) and ii) together with the monotonicity of Q state that the optimization 
problem (j33p is convex and admits a solution (for the sake of simplicity we also assumed differentiability of 
(p); whereas iii) guarantees that there exists a finite (large enough) r such that ^^^(n) ^C") is a Px NEP (cf. 
Lemma El]). The assumptions on the sequence {e^")}^]^ are also rather weak and require e^"^ to go to zero, 
but not too fast; which can be satisfied, e.g., by taking e^"^ = l/{l + na), with n = 0, 1, 2, . . . and a being any 
positive constant. The reason for this requirement is rather intuitive: if e^"^ became "too small, too soon", 
then the term ^||x.j — yj|p in the players' objective functions would dominate the term e^"^ (j), making the 
role of e^") (f> negligible; Algorithm [4] would become "almost" Algorithm |3] so that one could not guarantee to 
find nothing more than a point in SOL(Q, f). 

The implementation of Algorithm H] requires the ability of solving at each round n the P-f NEP G^ ^(„) . 
To this end, as already observed in the previous section, we can use several centralized methods. However, 
if we are interested in distributed solution schemes, Algorithm [T] applied to i,(n) ^(n) is the natural choice. 
Note that the convergence conditions of the algorithm applied to G^ ^(n) ^(n) ^s stated in Theorem [13] are 
always met, provided r is large enough; see, e.g., (j35p in Lemma [2T] 
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As in the previous section, note that, unless one has simple ways to compute the solutions of the games 
x(") 1 Algorithm [J] requires at each step the exact computation of the solution of the regularized games 
^Te(") x{") (inner loop), which in principle requires a conceptually infinite procedure. While in practice this 
might not be a problem, it still leaves open the question of whether, paralleling the results of the previous 
section, one can develop versions of Algorithm |4] where inexact solutions are used in the Step 2. The answer 
to this question is positive, but the corresponding theory is rather complex and, for sake of simplicity, we 
prefer to omit it here; the interest reader can work it out using results in [21) . 

5.4 A bird's-eye view 

In the previous three sections we proposed several distributed algorithms, which are applicable to different 
scenarios. It is useful to summarize the results obtained so far, showing that, in spite of apparent diversities, 
all the algorithms belong to a same family; Figure |3] provides the roadmap of the proposed distributed solution 
methods for real NE 



Real player-convex NEPs 

minimize f,(xi,x_i) 
subject to Xi G Qi 



Q 
is 



Q has a unique NE 



Q may have multiple NE 



Solve the Pt NEP Q using 
Algorithm 1 



Solve a sequence of P-r NEPs 
^?T.e(").x(") using Algorithms 1 




ST,e("),x(") 


: for each i = 1, . . . , /, 


mill /;(X;,x_ 


■ ) + I", ^x|"' 1 + e'"'0(xi.x_i) 


s.t. X,- e Si, 





^'r,el"),x("l • fol' each J = 1, . . . , J, 

T II ,,„,||2 

mill .fi(Xi, x_i) + - ||x, — Xj. II 
s.t. x,i e Si, 



Figure 3: The roadmap of the proposed distributed solution methods for real player-convex NEPs 

Conceptually, what we have proposed is indeed a unified algorithm, where the users can explicitly choose 
the degree of desired cooperation and signaling, converging to solutions having different performance, namely: 
i) any one NE, when there is no cooperation among users, and ii) the "best" NE (according to an outer merit 
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function (p), when there is some cooperation. The choice of one scheme in favor to the other as well as the 
merit function (j) will depend then on the trade-off between signaling and performance that the users are 
willing to exchange/achieve. The core of the proposed solution methods can be summarized in the following 
unified updating rule: at iteration n, the optimal strategy of each user i is 



x^ = argmm <; fi ( Xj, xli/ ) + vr^^ ' ( Xj, x^/ ) + ^ 



r 



x.-x^) 



(36) 



where the first term /j(xj,x^^) is the usual term in an iterative best-response algorithm, the second term 
7r^''"'*(xj, x^j^) (whose update is performed at iteration z/„,) can be interpreted as a nonlinear pricing in the 

is a (proximal) regularization. Observe 



r 

objective function of the users, and the third term — 



X,. - x(-) " 



that there are two iteration indexes: n is the main discrete-time unit, whereas Un is increased every few 
discrete-time units (e.g., if i/„ = [n/lOj , then i/„ is updated every 10 discrete-time units). 

The price function vr^''"-' can be interpreted as a measure of the "altruism/selfishness" of the users and 
represents the trade-off factor between signaling and performance. Indeed, we may have the following: 

• vr = (no cooperation): The users are not willing to cooperate; the best one can get is converge to 
any one solution of the game (i.e., with no control on the quality of the solution); this is guaranteed 
even in the presence of multiple equilibria if the NEP is monotone (Pt); 

• vr-''"^ 7^ (some cooperation): The users may exchange some signaling in the form of pricing through the 
function 7r^''"^(xj, xL7) = e('^")(/)(xi, xL7) and converge to the NE that minimizes the merit function 
0(x); convergence is guaranteed if the NEP is monotone. 

Remark 23 (Role of pricing) It is important to remark that the pricing term 'n'f'^^ does not need to be 
linear; moreover it has a specific and well understood role in the optimization of the system performance. 
This is a fundamental departure from current literature that uses linear pricing as a heuristic to improve 
the performance of a NE in power control games (see, e.g., |56) for scalar power control problems); in these 
works there is neither a proof of convergence of the modified game nor a theoretical explanation of the 
performance improvement due to pricing. As a direct product of our framework, we obtain instead a clear 
understanding of the meaning of the pricing; for example, a linear price in the form 

gK)^T withe(^") ^ 0, 

corresponds to the selection of a NE that minimizes the linear function ^r^Xj , resulting likely in better 
system performance. 



Finally, the proximal regularization in (|36p has the role to numerically "stabilize" the algorithm; the 
proximal gain r is the trade-off factor between the convergence stability and the convergence speed. The 
beneficial effects of a proximal regularization are well understood in the optimization literature; see, e.g., [53]. 



6 Variational Inequalities and Games in the Complex Domain 

All the results presented so far apply to real NEPs. However, in many applications, e.g., in digital com- 
munications, array processing, and signal processing, the variables involved in the optimization are complex 
numbers. For instance, in the MIMO problems introduced in Sec. [21 the optimization variables of the players 
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are complex matrices. For these applications, the reformulation of the problem into the real domain is awk- 
ward, and generally leads to final conditions that cannot be easily rewritten in terms of the original complex 
setup. It seems instead more convenient to work directly in the complex domain. This naturally calls for the 
definition of the VI problem in the domain of complex matrices, which is the aim of this section. The novel 
contribution in this direction is threefold, as listed next. 

Following the same approach as in the real case (cf. Sec. [3]), we introduce in Sec. 16.21 the minimum 
principle for constrained convex optimization problems in the domain of complex matrices, generalizing the 
already known complex gradient-vanishing conditions obtained in |40| for the unconstrained minimization of 
real- valued functions of complex variables. As an intermediate result, we also introduce a Taylor expansion of 
real-valued functions of complex matrices that is amenable to our MIMO applications; details can be found 
in Appendix [El The second important contribution is given in Sec. 16.31 where we provide the formulation 
of a VI problem in the domain of complex variables that can be seen as the natural generalization of the 
minimum principle introduced in Sec. 16. 2| where the (complex matrix) gradient is replaced by a (complex- 
valued matrix) function of matrix variables. We then introduce the concepts of monotonicity and P properties 
for the complex-valued matrix function of the VI and provide matrix conditions for this properties to hold. 
Interestingly, these conditions are the natural generalization of those obtained in Sec. 13. H provided that a 
new definition of Jacobian matrix tailored to complex-valued matrix functions is introduced. Finally, the 
third contribution is the connection between Vis and NEPs in the complex domain, as given in Sec. 16.51 
Our derivations are based on the use of M-matrix derivatives (also termed Wirtinger derivatives) that are 
preliminary introduced in Sec. 16.11 this section is motivated by the lack of a well-established notation and 
definitions for R-matrix derivatives of functions of complex (matrix) variables; two good tutorials on the 
subject (termed CM-calculus) are |40| 157) (note that |40| and |57| use a different formalism, we follow |40|). 

6.1 M-matrix derivatives 

In practical applications, we often deal with optimization of real-valued functions / : C 9 z i— )■ f(z) G M of a 
complex variable z that are not differentiable in C (termed also C-differentiable or holomorphic) |j However, 
the same univariate function / : C — t- M can also be viewed as a bivariate function of its real and imaginary 
components, i.e., /(z) = g{zR, zj), where : i— ?• R is a real-valued function of the real variables zr = Re(z) 
and zj = Im(z). This way, one may be able to replace the nonexistence of the C-derivative of / with the 
existence of the real partial derivatives of g{zji, zj), which is actually what one needs to compute a stationary 
point of the function. This motivates the introduction of the so-called R-derivative and conjugate R-derivative 
of / : C — >■ R at zq £ C, formally defined as 



df. . ^ 1 fdfiz) .df{z) 
(^o) = ^ J 



dz 2 \ dzR dzj 



, 5/ , 1 fdfiz) .df{z) 



(37) 



respectively, where j = Note that the derivatives above must be interpreted formally, because z 

and its conjugate z* in (|37p are treated as they were mutually independent; the derivatives and 
represent instead the true (non- formal) partial derivatives of / viewed as a bivariate function of zr and zj, 
i.e., / = f{zji,zj). When and exist (and are continuous), implying that ([37]) is well-defined, we say 



^It is a known fact that nonconstant real-valued functions (of complex variables) are not C-differentiable. 
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that / is R-difFerentiable (or continuously M-differentiable); similarly to the real case, when we say that a 
function / : /C — )• M is M-differentiable (or continuously M-differentiable) on the closed set /C, we mean that 
the function is so on an open set containing /C. 

The M-derivatives defined in (j37p for a real-valued function can be naturally extended to complex-valued 
functions of a complex argument, that is, / : C — t- C; formally we still have ()37p . but now f{z) = f{zji,zi) = 
fRizR,zi) + j ■ fi{zji,zj), with / : M2 C and /rJi : M^ h-^ M, and by df/dzR we mean df/dzR = 
dh/dzR+j ■ dfi/dzR (similarly for df/dzi). 

When / is a (complex-valued) scalar function of complex matrices, that is / : C"^'" — t- C, we have n ■ m 
component-wise M-derivatives ^^^-^ and n- m conjugate M-derivatives q^^^*) ■ '^^^ question naturally arises 
how to order these n ■ m complex derivatives. Obviously this can be done in many ways, corresponding to 
the pattern chosen to arrange the n ■ m elements. It is worthwhile noticing that, even though they all contain 
the same n-m derivatives, not all definitions have the same properties; for instance for some of them a useful 
chain rule does not exist. Next, we introduce two definitions, both useful for our derivations and widely used 
in the literature |40| : in the former definition, the n • m (conjugate) M-derivatives are displayed in the same 
order as (Z)jj and (Z*)jj appear in Z and Z*, whereas in the latter we arrange all the elements in a row 



vector. Given / : C"^™ 
as 



the (matrix) gradient and co(njugate)-gradient of / at Zq S 



are defined 



Vz/(Zo) 
Vz*/(Zo) 



df{Z) 



dZ 

A df{Z) 



Z=Zo 



dZ* 



with 



with 



z=z 







p/1 




df 


dZ_ 




diZ)ij' 


\df 




df 


dZ* 




d{z*)ij 



Vi = 1 , . . . , n and j = 1, . . 



yi = 1, 



and j = 1, 



, m 



(38) 



where ^^^^ and are the M-derivative and conjugate M-derivative of the complex- valued function / 



9(Z),, 
w.r.t. (Z)-- and 



Z*)--, respectively. Note that Vz/(Zo) and Vz*/(Zo) are matrices having the same size 



of Z. Alternatively, one can arrange the elements 
L'z*/(Z) at Zo G V 



df 
9(Z), 



and 



df 
9{Z*),, 



in a row vector, and define Dzf{Z) and 



as 



DzfiZo 



dfiZ) 



dvec{Zy 



= vec(Vz/(Zo))^ andZ)z*/(Zo)^ ^ ^ 



dfiZ) 



Z=Zn 



vec(Vz./(Zo))^ , 
(39) 



where vec (A)"^ stands for (vec(A))'^. For (complex-valued) matrix functions of complex matrices, F*^ : 
(^nxm _j, (j^pxg^ arrange the pq ■ nm (conjugate) M-derivatives in the following pq x nm matrices 



L)zF'^(Zo) 
Dz.FC(Zo) 



dvec (F'^(Z)) 



9vec(Z) 
9vec (F'^(Z)) 



with 



z=z 







9vec(Z 



with 



Z=Zo 



5vec (F*^)' 
dvec{zf 

dvec (F*^) 
dvec{Z*f 



\/i = 1, . . . ,pq and j = 1, . . . , nm, 

, yi = 1, . . . ,pq and j = 1, . . . , nm. 

(40) 

The pq x nm matrices DzF"^ and Dz*F'^ are called Jacobian and conjugate Jacobian of F*^. Note that when 
F^ is a scalar function of Z, i.e., F^(Z) = /(Z) with / : C"""™ C, definitions (gO]) reduce to ([M])- Practical 



ghc(F^)], 
9[vec(Z)]^. 

d [vec(FC)]^ 
" d [vec(Z*)]. 
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rules to compute M-derivatives and conjugate M-derivatives introduced above for specific functions of vector 
and matrix variables can be found in |40|. 



6.2 The minimum principle 

Let C"^™ be the space of complex n x m matrices, and let /C C (C"-x"i be a closed and convex set. We 
consider the optimization problem 

minimize f(Z) 

(41) 

subject to Z € /C, 

where / : /C — t- M is a real- valued convex and continuously M-differentiable function on /C. At the basis of the 
minimum principle there is the first-order Taylor expansion of / at Zq E /C as proved in Appendix [eI^ 

/(Zo + AZ) - /(Zo) ~2Re(tr ((Vz/(Zo)f Az) ) 
= 2Re(tr(AZ^ (Vz./(Zo)))) 

^2(AZ, Vz./(Zo)) (42) 

where we used (Vz/)* = Vz*/ since / is real [see (|37p ]. and in the last equality we introduced the inner 
product (•, •) : C"^" x C"^" ^ M, defined as 

(A, B) ^Re(tr(A^B)). (43) 

Note that the norm induced by the inner product (•, •) is the Frobenius norm, i.e., (A, A) = Tr(A^A) = 
II A||^ . Using (j42p we can now introduce the minimum principle as given next. 

Lemma 24 Given the convex optimization problem in the setting above, 'K & IC is an optimal solution 
of (^^d. only ifX satisfies (Z - X, Vz*/(X)) > for all Z G /C. 

Proof. See Appendix [E) ■ 

It is interesting to observe that if the optimal solution X is in the interior of /C [e.g., the optimiza- 
tion problem (j4ip is unconstrained, implying /C = C"'^™'), then the above optimality conditions reduce to 
Vz*/(X) = 0, or equivalently Vz/(X) = 0, which are the well-established complex gradient- vanishing con- 
ditions obtained in |40| for the unconstrained minimization of real functions of complex variables. Here, we 
generalize that result to the case of constrained minimization of real-valued functions of complex matrices. 
In fact. Lemma [23] provides a platform for studying constrained optimization problems within the complex 
variables domain, taking advantage of the multivariate CM-calculus (and in particular complex matrix deriva- 
tive rules |40| ) as well as the structure of the problem under consideration. We conclude this section with an 
example of application of the minimum principle, which is instrumental for the analysis in Sec. [71 



^The proposed Taylor expansion can be rewritten as the one reported in [401 Ch. 5], for which however there is no proof. 
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Example 25 (An application of the minimum principle) Consider the fohowing single-user rate max- 
imization problem 

maximize /(Z) = log det (R^ + HZH^) 

z (44) 
subject to Z G /C, 

where R„ G ^jmxm -g ^ positive definite matrix, H G C"^^", and /C is any convex and compact subset of 
the nxn complex positive semidefinite matrices (assumed to be nonempty). Note that /(Z) is a concave 
(real- valued) function on the feasible set /C; however, / is not real if defined on C"^". Since we are interested 
in minimizing / over /C, in order to apply the minimum principle, one approach we can follow is to consider 
without loss of generality the modified function / : /C — )■ M, defined as /(Z) = 2 Re (/(Z)), where /C C C'^x'^ 
is any open set over which /(Z) is well defined (it is sufficient that det (R„ -|- HZH^) 7^ 0); indeed, 
coincides with / over /C, but it is real everywhere (in its domain). Moreover, / is M-differentiable on /C. 
The conjugate (matrix) M-derivative of / at Zq G /C is (see Appendix |F|) 

Vz. /(Zo) = H^(Rn + HZ^H^)-iH. (45) 

Introducing the complex- value matrix function G : /C — )• C"^", defined as G = G(Z) = — Vz*/(Z) = 
-H^(R„ + HZH^y^U (note that Z G /C and thus Z = Z^), and invoking Lemma [231 the optimization 
problem (j44p is then equivalent to the minimum principle in the complex domain: find a Z G /C such that 
(Y - Z, G(Z)) > for aU Y G /C. □ 

6.3 The VI problem in the complex domain 

With the developments of the previous section on hand, we can now introduce the definition of the VI problem 
in the domain of complex matrices, termed the complex VI problem. Similarly to what we have done for the 
real case (cf. Sec. [3|), we can think of the VI problem as the generalization of the minimum principle (cf. 
Lemma [24|) . where the co-gradient Vz*/ is replaced with a complex- valued matrix mapping. The formal 
definition is given next. 

Definition 26 Given a convex and closed set K, C C"x™ and a complex-valued matrix function F'^(Z) : /C 3 
Z —7- C"^"^, the complex VI problem, denoted by VI (/C,F^), consists in finding a point Z G /C such that 
(Y - Z, F^(Z)) > for all Y G /C. The solution set of the VI (/C, F^) is denoted by SOL (/C, F^) . 

When the the feasible set /C has a Cartesian structure, i.e., /C = 111=1^* with each /Cj C c^-ixm-i^ 
write F^(Z) = (Ff(Z))J^^ and Z = {Zi)l^^, with Ff(Z) : /C ^ C"'Xm, r^. ^ ^n,xm,^ ^^^^i a 



*The introduction of the auxiliary function / might appear an unnecessary complication, which needs clarification. The 
original function / is defined over a (sub)set of positive semidefinite matrices. The theory of matrix derivatives introduced in 
this paper cannot be applied directly to functions of matrices having a structure. The function / allows us to overcome this 
issue, because it is defined over an open set of unpatterned matrices (over which we can use the developed matrix differentiation 
theory) while being equal to / over the set of interest. A different approach would be working directly with the original function 
/ and using the so-called complex (patterned) generalized derivatives |40] . However, state-of-the-art results show that it is 
not possible to find generalized derivatives w.r.t. an arbitrary set of complex-valued patterned matrices, which strongly limits 
the applicability of this methodology in practice. Using instead the approach explained in this example, the proposed matrix 
differentiation machinery can be applied also to optimization problems involving patterned matrices. 
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case, with a slight abuse of notation, we will still use for the partitioned VI (/C, F''-') the compact notation 
(Y - Z,F^(Z)) > 0, by meaning Yli=i (Yi - Zi,F'^i(Z)) > 0. Moreover, the definitions of DzF^{Z) and 
Dz*F*^(Z) as given in (|40p depend in principle on the ordering according to which the components of F'^(Z) 
and Z are grouped in the vec operator. For our purposes, the following ordering is the most convenient, 
which is tacitly assumed throughout the paper: vec ((F^ (Z))[^-^) = [vec(Ff (Z))-^, . . . , vec(Fj (Z))^] and 
vec ((ZOLi) = [vec(Zi)^, . . . , vec(Z,) 



T" 



T 



6.4 Monotonicity and P properties of VI (/C, F*^) 

We can now readily extend the definitions of monotonicity and P property as introduced in Proposition [16] 
to the complex- value matrix mapping F^; the aforementioned definitions are in fact formally the same, with 
the only difference that the scalar product and the Euclidean norm are replaced with the inner product (•, •) 
defined in (j43p and the Frobenius norm, respectively. The non-trivial task is instead to provide conditions 
easy to check that guarantee these properties. These conditions are indeed instrumental to study convergence 
of algorithms for complex NEPs. The interesting result we prove next is that we can obtain necessary and 
sufficient conditions for a continuously (]R-)differentiable F'^(Z) to be a monotone function or (sufficient 
conditions to be) a P function on /C that are formally equivalent to those obtained for real-valued vector 
functions F(x) [cf. (f20|) and Proposition [5] , provided that we introduce a novel definition of Jacobian matrix 
suitable for complex-valued functions of complex variables; such a Jacobian will contain both M-derivatives 
and conjugate M-derivatives of F'^(Z). 

Given the complex VI (/C, F"^), suppose that F'^(Z) is a continuously (]R-)differentiable matrix function 
on /C. Then, the nm x nm Jacobian matrices i:>zF^(Z) and Dz'F^{Z) in (gO]) are well-defined at Z G /C. 
Let us introduce the 2nm x 2nm matrix JF'^(Z), defined as 

Z)zFC(Z) Dz*FC(Z) 
Dz (FC(Z)*) Dz* (FC(Z)*) 

which we call "augmented Jacobian" for obvious reasons. For notational simplicity, in the sequel we will write 
DzF'^(Z)* and L»z*F^(Z)* for Dz (F^(Z)*) andL>z* (F'^(Z)*), respectively. Note that the following rela- 
tionships hold between the blocks of JF'^(Z): (DzF^(Z))* = Dz*F^(Z)* and (Z)z-F'^(Z))* = DzF^{Z)*. 
Finally, under the assumption that F'^(Z) and }C have a partitioned structure and F'^(Z) has bounded 
(M)-derivatives on /C,, let us introduce the "condensed" I x I matrix Tpc given by 



JFC(Z) ^ - 



(46) 



t""^^]-- - ^ otLrwise, ^''^ 



with 

k: 



^ inf Aieast Af J,Ff(Z)Ai and " Af J,Ff (Z) A,- 



F ' 



Z6/C V ' " / z&K 

where JjF^(Z) and JjFp(Z) represent the augmented Jacobian matrices of F^(Z) as defined in (j46p . whose 
M-derivatives are taken with respect to the matrix variables Zj and Zj (and their conjugates), respectively; 
Aj G (£2nimix2nimi ^^.^ nonsingular arbitrary matrices; and ||A||^ denotes the Frobenius norm of A. Matrices 
JF*^(Z) and Tpc play for complex Vis the same role as JF and Tp introduced in Sec. 13.11 for real Vis. 
Monotonicity and P properties of VI (/C, F"^) are indeed relayed on those matrices, as shown next. 
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Before stating the main results (Propositions [28] and I29p . we need to introduce a novel relaxed definition 
of (uniformly) positive (semi-)definiteness for matrices in the form (|^6|) . which takes explicitly into account 
the special structure of those matrices. Instead of checking the sign of the quadratic form y^JF'^(Z)y for 
all possible vectors y G C^*^*", it turns out that we can restrict the check to structured vectors in the form 
y = [yi,yi] for all yi G C"*", which is actually the size of the vector space where Z lies (see Propositions [28l 
andl29lbelow). This motivates the following definition of "augmented" (uniformly) positive (semi-)definiteness 
for matrices in the form of (|46p . 

Definition 27 The augmented Jacobian JF'^(Z) is said to be: 

i) augmented positive semidefinite on /C if for all Y G C^x™ and Z G /C, 

vec([Y,Y*])^JF'^(Z)vec([Y, Y*]) > 0; (49) 

ii) augmented positive definite on fC if for a// 7^ Y G C"^™ and 7^ ^ K,, the inequality in ( |-^9[ ) is strict; 

iii) uniformly augmented positive definite on fC with constant c > if for all Y G C"^™ and Z G /C, there 

exists a positive constant c > such that 

vec([Y, Y*])^JF^(Z) vec([Y, Y*]) > c ||Y||^ . (50) 

For i), ii), and iii) we will write JF^(Z) ^ 0, JF^(Z) y 0, and JF'^(Z) - cl ^ 0, respectively. 

Note that JF'^(Z) is not Hermitian; which implies that vec(W)-'^ JF'^(Z)vec(W) is generally not a real 
number for arbitrary W G C"^^*". However, because of the structure of JF^(Z) and vec([Y,Y*]), with 
Y G C"'^"^, the quadratic form vec([Y, Y*])^JF^(Z) vec([Y, Y*]) introduced in the proposition is always 
real. Note also that if JF'^(Z) is positive semidefinite, positive definite, or uniformly positive definite on /C 
(and thus Hermitian), then it is also augmented positive semidefinite, positive definite, or uniformly positive 
definite, respectively; but the converse in general is not true (because JF'^(Z) is not Hermitian). Using this 
new relaxed concept of positive definiteness, we can now establish the connection with the monotonicity 
properties of F^. 

Proposition 28 Let F^ : /C— t- C"^"' be (R-) continuously differentiable on the convex set /C. Suppose that K, 
has nonempty interior. The following statements hold: 

(a) F*^ is monotone on K. if and only if JF'^(Z) is augmented positive semidefinite on /C; 

(b) // JF'^(Z) is augmented positive definite on K., then F*^ is strictly monotone on K,; 

(c) F*^ is strongly monotone on K, with constant Cgm > if and only JF'^(Z) is uniformly augmented positive 

definite on K, with constant Csm/'^- 

If we assume a Cartesian product structure, i.e. F^ = {Ff)l^^ andKL = 11^=1 '^^'^ bounded (M.) -derivatives 
o/F*^ on IC, then: 
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(d) 7/Tpc is positive semidefinite/PQ-matrix, then F*^ is a monotone/Po function on K,; 

(e) //Tpc is a P-matrix, then is a uniformly P-function on /C. 

Proof. See Appendix [Gl ■ 

The above proposition is the generahzation of (j20p and Proposition [5]to complex Vis. Note that, if the set 
JC has empty interior, necessary conditions in (a) and (c) generally do not hold, whereas sufficient conditions 
in (a)-(c) may be too restrictive. Since some of the optimization problems of our interest have feasible sets 
that fall into this class [e.g., think of the set of Hermitian matrices], it is worth extending Proposition [28l to 
sets with empty interior. The next result is valid for arbitrary (nonempty) convex sets. 

Proposition 29 Consider the setting of Proposition \28\ but with IC being any nonempty convex subset of 
(£^nxm^ Lei Sfc be the subspace that is parallel to the affine hull of The following statements hold: 

(a) F*^ is monotone on JC if and only if for all Y € £ny.m ^^^/j ^/^^^ y g Sjc and Z ^ IC, it holds 

vec([Y,Y*])-^ JF^(Z)vec([Y,Y*]) > 

(b) If for allO^Y £ C"^™ such thatY £ Sk: andZ £ IC, it /loZds vec ([Y, Y*])^ JF^(Z) vec ([Y, Y*]) > 0, 

then F*^ is strictly monotone on IC; 

(c) F*^ is strongly monotone on IC with constant Cgm > ^ if oind only for all Y G £;"x»" such that Y G Stc 

andZeIC, it holds wec{\Y ,Y*])^ 3¥^{Z)wec{\Y ,Y*]) >{c^^/2) \\Yfp. 

If we assume a Cartesian product structure, i.e. F*^ = (Fp)^^-^ andIC = Y[i=i^i! '^^'^ bounded [M.) -derivatives 
o/F'^(Z) on IC, then statements (d) and (e) of PropositionWR hold. 

Proof. See Appendix [Gl ■ 

A set IC of special interest for our applications (cf. Sec. [2]) is the set of complex n x n positive semidefinite 
matrices (and thus Hermitian). This set has empty interior, implying that one needs to use Proposition [29l 
It is not difficult to see that the affine hull of such a /C is the set of Hermitian matrices, which is already a 
subspace. Therefore, when Proposition 1291 applies to such a /C, the matrices Y are restricted to the set of 
Hermitian matrices. It is worth observing that, when IC has nonempty interior. Proposition [29] reduces to 
Proposition [28] indeed, we have Aff(/C) = C"^™, and thus Sk = C"^*". 

Using Proposition [28] (or Proposition [29]) and building on the structure of JF'^(Z) one can obtain sufficient 
conditions for JF^(Z) to be augmented positive (semi-)definite or uniformly positive semidefinite, similarly 
to what we have done in Sec. 13.11 for real valued vector functions F(x); one can then extend the solution 
analysis and methods developed for the VI(Q, F) to the complex VI (/C, F*^) ; because of the space limitation, 
we leave these tasks to the reader. In Sec. [7] we will show an instance of these conditions when specialized 
to the MIMO games along with their physical interpretations. 

We conclude this section by applying Proposition [28] (or Proposition [29]) to the conjugate gradient of 
real- valued functions of complex variables [cf. (|4ip ]. The result is a set of novel necessary and sufficient 



^We recall that, given a subset K, of C"^"", the affine hull of /C, denoted by Aff(A3), is the set of all affine combinations of 
elements in K., that is A&{]q = |y £ C"^™ : Y = ^Li fc > 0, X, G /C, a, £ R, ELi = l}- 
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conditions for a (continuously M-differentiable) real-valued function of complex variables to be (strictly) 
convex or strongly convex, in terms of R-derivatives. This provides an easy way to check convexity directly 
in the complex domain. In order to apply Propositions [28] or [29l we need the following intermediate result, 
which can be proved using the Taylor expansion (|42p and following the same approach used to prove similar 
results for real- valued functions of real variables. Given a continuously M-differentiable real- valued function 
/ : C^x"^ — ]R , / is convex, strictly convex, or strongly convex on /C if and only if its conjugate gradient 
Vz*/(Z) is monotone, strictly monotone, or strongly monotone on /C, respectively. Using Proposition [28] 
(or Proposition [29l) . the convexity properties of /(Z) can be then restated in terms of properties of the 
augmented Jacobian matrix JF'^(Z) of F^(Z), with F*^(Z) = Vz*/(Z), which we term augmented Hessian 
of /(Z), nzz*f{Z), given by [cf. m]- 



Hzz*/(Z) ^ ^ 



I?z(Vz*/(Z)) Dz*(Vz./(Z)) 
DziiVz* f{Z)r) Dz^iiVz' f{Z)r 



1 



V|z./(Z) V|.z,/(Z) 
/(Z) V|.z/(Z) 



^ZZ 



(51) 



Note that [cf. 1^] V|z./(Z) = (V|,z/(Z))* and V|,z,/(Z) = (V|z/(Z))*. It follows from Proposition 
[28] applied to F^(Z) = Vz*/(Z) that ^zz*/(Z) plays the role of the classical Hessian matrix of /: let 
/C C C"x™ be any convex set with nonempty interior, then 



/(Z) is convex on /C 
/(Z) is strictly convex /C 
/(Z) is strongly convex on /C 



A 



^ ^zz-/(Z) ^0, VZg/C; 



A 



nzz*f{Z)yO, VZG/C; 



(52) 



A 



'Wzz*/(Z) — CsmI ^ 0, VZ G /Cand some Cgm > 0. 



If the set /C has empty interior, conditions (I52p are replaced by those in Proposition [29] applied to 
JF'^(Z) = ^zz*/(Z); we leave this easy task to the reader. Note that our conditions in (j52p (and Proposition 
[29]) generalize those obtained in |58| Prop. 1.2.6 and Exercise 1.8] for real-valued functions of real variables, 
and sufficient conditions reported in |40| Ch. 5.3.1] for the convexity (only) of real- valued functions of complex 
variables defined over C"^™ (i.e., /C = C"^™). 

Example 1251 Revisited. Going back to the optimization problem (I44p . we can recover the well-known 
concavity property of /(Z) on the compact and convex set /C by a direct application of Proposition [29] The 
expression of the augmented Hessian of /(Z) will be also used in Sec. 17.21 for studying MIMO games. 

Let JC be any open set containing /C over which /(Z) is well de&ied, and let / : /C — )• M be /(Z) = 
2 Re (/(Z)). Since / = / on /C, concavity of / on /C follows from that of / on /C. Since /C has empty interior, 
one needs to use Proposition [29] Observing that, for the specific set /C under consideration, the set Sjc in 
Proposition [29] is ^yc = {X G C"^" : X = X^}, it is sufficient to show that 

H 



vec(Y) 
vec(Y*) 



vec(Y) 
vecfY*) 



> 0, VZG/CandVY = Y 



nzz'fiZ) 

where 'Hzz*/(Zac) is the augmented Hessian of /(Z). In Appendix [F] we show that 

nzz*f{z) = 



H 



(53) 



[G(Z)^' 





G*(Z) 



K 2 2 



[G(Z)^®G(Z)] K„2 




(54) 
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where G(Z) = (R„ + HZ^H^) ^ H and K„2„2 is an x r? permutation matrix such that vec(Z'^) 
K„2„2vec(Z) (also termed commutation matrix [401 Def. 2.9]). Using (|54p . condition (|53p becomes 



< 



vec(Y) 


11 p 


vec(Y*) 




vec(Y) 


H r 


vec(Y*) 




vec(Y) 


H r 


vec(Y*) 







G(Z)-f^^G(Z)* 

G(Z)^»G(Z) 


G(Z)^»G(Z) 




G(Z)^®G(Z) 








G(Z)-f^ ® G(Z)* 







G(Z)^ 8)G(Z)* 



K„2„2 vec(Y) 
K„2„2 vec(Y*) 

vec(Y^) 
vec(Y^) 

vec(Y) 
vcc(Y*) 



VZ G /CandVY = Y 



(55) 



where in the second equahty we used the property K„2„2vec(Z) = vec(Z'^), wheres in the last equality 
the condition Y = Y^. It turns our that (|55|) is satisfied if G(Z)-^ ® G(Z) is positive semidefinite for 
all Z G /C. Since G(Z)-^ ® G(Z) is Hermitian on /C, it is sufficient to check that the minimum eigenvalue 
of G(Z)^ ® G(Z), denotes by Amin (G(Z)^ ® G(Z)), is nonnegative for all Z on /C. The result follows 
from Amin (G(Z)^ ® G(Z)) = A^in (G(Z)^) • Amin (G(Z)) = Amin (G(Z))2 > for all Z G /C, which proves 
concavity of /(Z) on /C and thus of /(Z) on /C. 

It is worth observing that while — 'Hzz*/(Z) satisfies ([53|) [Proposition [29lfa)]. it is not positive (semi-) 
definite, showing that the latter condition may be too restrictive for checking the convexity of a (real-valued) 
function of complex variables. This strengths the importance of the proposed new concept of augmented 
positive (semi-)definiteness (Definition and the role of Propositions [251 and [201 □ 



6.5 NEPs in the complex domain 

We can now establish the formal connection between complex NEPs and complex Vis. Let t/c — (^i f ) be 
a complex NEP where each player controls a complex matrix Zj G ([^nixm.i ^j^^^^ must belong to the player's 
feasible set /Cj C the cost function of each player is denoted by /j : /C — )• M; and the joint strategy set 

of the game is /C = Hi ^i- We also write Z = (Zj)^,!, Z_j = (Zi, . . . , Zj_i, Zj+i, . . . Zq), and /C_j = Ylj^iJCj. 
The NEP problem Qc consists then, for each player i = 1, . . . , Q, in solving the following convex optimization 
problem: given Z_j G IC-i, 

minimize /i(Zj,Z_j) 

(56) 

subject to Zj G /Cj. 

Building on Lemma [2H it is not difficult to prove the following. 

Proposition 30 Given the complex NEP Q^. = {IC,i) , suppose that for each player i the following hold: 

i) the (nonempty) strategy set /Cj is closed and convex; 

ii) the payoff function /j(Zj, Z_j) is convex and continuously M-differentiable in Zj for every fixed Z_j. 
Then, the complex NEP is equivalent to the complex Vl{IC,G'^), where G'^(Z) = (Vz*/j(Z))^_i^. 

Exploring the above connection between Qc and the VI(/C, G^) as well as the monotonicity/P properties of 
VI(/C, G^) (cf. Proposition 1281 and Proposition I29p. one can then extend the solution analysis and algorithms 
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developed for the real NEP Q m Sec. [3l U and [5] to the complex NEP [VI(/C, G"')]; because of space 
limitations, we leave this task to the reader. In the next section we will show how to explore the connection 
between Qc and the VI(/C, G^) along with the developed theory of complex Vis to study and solve the MIMO 
games introduced in Sec. 12.21 



7 Noncooperative Games Over Interference Channels Revisited 

In this section, we focus on the application of the general theory developed in the previous sections to some 
concrete examples of practical interest. In particular, we show how the real/complex NEPs introduced in 
Sec. [2] can be naturally casted in the proposed framework and thus efficiently solved. The main result 
in the SISO case is a novel iterative water-filling like algorithm where the users can choose the degree of 
desired cooperation via local pricing, converging to solutions having different performance/signaling trade- 
off; we also prove that the best-response of each player has a multi-level water-filling-like expression and 
provide an efficient algorithm for its computation. We then extend our analysis to MIMO games and obtain 
similar results. Numerical experiments show the superiority of our novel algorithms with respect to plain 
noncooperative solutions as well as very good performance with respect to centralized schemes, in favor of 
very limited signaling among the users. 

7.1 The SISO case 

We study here the game ^siso — ^-^siso^ (^2)^=1) ii^troduced in (iTl) . Tlie VI function associated witli ^siso 
is G(p) = (Gi(p))f^i : V^^^° R^^ , where each Gi(p) is defined as 

Note that in this section, due to the nature of the problems at hand, we called the VI mapping G instead of 
F used previously, and the VI variables p instead of x as used previously. 

According to Proposition [5l the monotonicity/P properties of G(p) are related to the matrices JGiow 
and Tg defined in ()2ip and (j22p . We recall that the matrices B and Cj's in the definition of JGiow and Tq 
represent a degree of freedom that one can use; in this case it is convenient to make the following choices. 
Let us rearrange the components of p by subcarriers, meaning that the vector p = (pj)f=i is permuted into 
p = (p(fc))^x, with p{k) = {pi{k))l^i; it is not difficult to see that p can be written as p = Pp, where P is 
a permutation matrix defined as 



1, if i = [(imod I)-l]N + mod(I • N) 
0, otherwise . 



Using this new order of the variables, matrix JG becomes P-'^JGP; JGiow is then obtained from P-^JGP ac- 
cording to (j2ip . where B = Diag |(B(A;))^^| is a block diagonal matrix, with each block B(A;) G R^^^ being a 

positive diagonal matrix with the i-th entry equal to [B{k)]-- = a^{k)/\Hii{k)\^+ Ejd^ij (^)P/l^^^(^)P)P™''''(^)• 
Matrix Tg comes directly from the original JG by choosing each Cj G M^^^ as a diagonal matrix, defined 
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as C,- = Dia; 



g I ((^?(^) + J2j \Hij{k)\^pf^-{k)'^ J. The explicit expressions of JGiow and 

are the fohowing: JGiow — Diag |(JGiow(^))fcLi| ^ is ablock diagonal matrix, whose k-th diagonal 



block JGiow(fc) G K'''' is 



A 



^1, if i = i 

[ JGlow(^)]ii 



and Tg G M is given by 



H.m' ■ ...^ . k = l,...,N, (58) 

^-^•.nn,,(fe), if^/,, 



if i = j 



\H^m' ..1 .... (^9) 

l<fc<7V 1 |i; 

with 



— max < — 7T ■ innr,,(A;) > if i 7^ 7, 



mmij{k) = 2— . (60) 



^j{k)+E r\Hjr{k )\^P?^''{k) 

Using the above definitions along with Proposition [H Theorem [H and Corollary [101 the main properties of Q 
are then given in the following proposition. 



Proposition 31 Given the real NEP Gsiso = (V^'^^°, (ri)[^^), the following hold. 

(a) Gsiso is equivalent to the Y1{V ^G), which has a nonempty and compact solution set; 

(b) Suppose that JG^ow is positive semidefinite (positive definite). Then G(p) is monotone (strongly mono- 

tone) on J> ; therefore Gsiso is a monotone NEP; 

(c) Suppose that Tg is a P-matrix (positive definite matrix). Then G(p) is a uniformly P-function (strongly 

monotone function) on 'p^'^^°; therefore Gsiso is a P-f NEP and has a unique NE. 

Using Proposition [TOl we obtain the following corollary providing sufficient conditions for the matrices in 
Proposition 1311 to be P or positive definite. 

Corollary 32 The matrix Tg in (|59p is a P-matrix (or a positive definite matrix) if one (or both) of the 
following two sets of conditions are satisfied: for some w = (wj)^^^ > 0, 



Low received MUI: — ' / ™ax_ \ , ■,,9 ' innrj,(/c) ^ < 1, Vi = I,-- - ,1, 



1 V- i \H^i(k)\^ 
— > Wj max < ^ • innr,o(K) > < 1 

^ 'l<k<N\\H,,{k)\' 'j 

Low generated MUI: — \^ Wi max < 

Wi l<k<N 



J^^^^^-innr,,(A;)l <1, Vj = l,---,I 



(61) 



Similar sufficient conditions can be obtained for JGiow (^) to be positive semidefinite. 

These conditions have an interesting physical interpretation: the uniqueness of the NE is ensured if the 
interference among the SUs is sufficiently small; this is clearly shown by Corollary [32j Specifically, the 
first condition in ()6ip can be interpreted as a constraint on the maximum amount of interference that each 
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receiver can tolerate, whereas the second condition introduces an upper bound on the maximum level of 
interference that each transmitter is allowed to generate. We will show shortly that these conditions play a 
role also in the convergence of the proposed distributed iterative algorithms. Moreover, depending on the 
level of interference in the network, the NEP ^siso is a P-f or monotone NEP, implying different properties 
and solution schemes of the game, as described next; we classify these two scenarios as low-interference and 
medium/high-interference regime, respectively. 

The case of Pf NEP (low- interference regime). When the matrix Tg is a P matrix (or positive definite), 
the NEP Gs iso is a Px NEP [Proposition [ST] (c)]. Invoking Theorem [ISJ the unique NE of the game can be 
computed with convergence guarantee using Algorithm [1] on t/siso as stated in the next theorem. 

Theorem 33 Suppose that Gsiso is a P-f real NEP. Then, any sequence {p^"^}$^o generated by Algorithm 
[7] applied to Gsiso converges to the unique NE of the NEP, for any given updating schedule of the players 
satisfying assumptions Al)-A3). 

When implementing Algorithm [H each user needs to compute his best-response solution, given the in- 
terference generated by the others. In Sec. 17.1.11 (cf. Lemma [35|) . we prove that the best-response for the 
game ^siso has a multi-level waterfilling-like expression, implying that each user can compute his optimal 
solution locally and very efficiently (he only needs to measure the overall MUI experienced at his receiver 
and "waterfiH" over it). Therefore, Algorithm [1] results to be totally distributed and computationally efficient, 
which makes it appealing for practical implementation in CR scenarios. 

The case of monotone NEP (medium/high-interference regime). When JGiow is positive semidefinite, 
the NEP Gs iso is a monotone NEP, having in general multiple solutions. In such a case, to compute a 
solution of Gs±so with convergence guarantee, there are two available options, namely: PDAs (either in their 
exact or inexact form) and PTRA. The former are the only feasible choice when the SUs are not willing to 
cooperate; whereas the latter requires some (albeit very limited) signaling among the SUs in favor of better 
performance (the algorithm converges to the "best" available solution). To the best of our knowledge, the 
above algorithms are in the signal processing and communication literature the first example of distributed 
power control schemes that converge even in the presence of multiple Nash equilibria. Note that, in all 
the aforementioned algorithms, the best-response of the SUs can be efficiently computed via a multi-level 
waterfilling expression; we address this issue in Sec. 17. 1.11 below. We provide next an instance of the PTRA 
along with its convergence conditions; PDAs are obtained as special cases of the PTRA, we thus omit their 
description. 

Equilibrium Selection via Proximal- Tikhonov Regularization Algorithm. According to the developments of 
Sec. 15.31 the first step is to choose a merit function that quantifies the quality of a NE of t/siso- Different 
heuristics can be used; as an example, here we focus on the following merit function: given the vector 
w ^ {wi)l^^ > 0, let 

/ TV 

1=1 jj^i k=l 
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This choice is motivated by the intuition that among all the Nash equilibria of ^, a good candidate is the 
one that minimizes the overall interference among the users, measured by the function (j){p), likely resulting 
in an "higher" sum-rate X]f=i^i(p)- The NE selection problem based on the merit function (p can be then 
formulated as: 

minimize i?^(p) 

P (63) 
subject to p G SOL(P^i^°, r). 

Problem (I63p is an instance of (j33p : we can then solve it by applying Algorithm [J] introduced in Sec. 15.31 

which corresponds to solving a sequence of perturbed P-f NEPs denoted by ^(„) p = (|7?=i^°^ (^_7.^(p)_l_ 

I 

p, and > 0, 

maximize nipi, p_i) - e^") 'jfpi - - ||pi - Pi||^ f64) 

AT A / \^ 

where 7 = (7j)f=i, with each 7^ = ( X^j^j ] • A partially asynchronous version of Algorithm |4] 

applied to (j63p is described in Algorithm [5] below, and its convergence conditions are given in Theorem [3~ 
which is a direct application of results in Theorem [13] and Theorem 



gC") 'yjp- + — IIp^ — pj||^^ y wherein each player i solves the following optimization problem: given p_ 



Algorithm 5: NE selection for the real NEP ^sis 



(Data) : {e^} ^ and r > 0. 



(S.O) 
(S.l) 
(S.2) 



Choose any p^'^^ G -psiso center p > of the regularization; set e = e^^^^ and n = 0. 

If p(") satisfies a suitable termination criterion, STOP. 
For each i = 1, . . . , I, compute p -""^^^ as 



p* G argmax In (pi,p^i '^"-"M - ejjpi - ^ ||pi - Pif | , if n G 7^ 

_ ^ ■T> siso ^ ^ / Z .* 

Pit / j 

pj""* , otherwise 



(65) 



(5.3) : If p("+^) is a NE of ^(„) p, then update e and the center p: 

£ = and p, = p["+^) Vi = l,...,/; (66) 

(5 . 4) : n n + 1 and return to (S.l). 



Theorem 34 Given the optimization problem (j63p . suppose that: i) Qsiso is a monotone NEP; ii) {e^"")} is 
such that e^"^ — )■ and X^^g ~ '^'^'^ ^^^v' ''" ^■^ chosen so that Tg + rl is a P matrix. Then, the 
sequence {p''"^}^o generated by Algorithmic has a limit point and every such limit point is a solution of 

A sufficient condition for matrix Tq + r I in Theorem [34] to be P is 



T > 
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Algorithm [5] shows that, in the presence of multiple equilibria, one can still have converge even when 
best-response based schemes (cf. Algorithm [1]) fail, provided that the SUs play a "sequence" of games rather 
than a one-shot game; moreover, to reach the NE that minimizes the overall MUI among the users, the 
players' objective functions need to be modified in order to contain an additional term— the linear term 
e^^^'y^'pj— whose task is to "measure" on the way the quality of the solution that the algorithm is going to 
reach. Such a term has also a physical interpretation: it represents a punishment imposed to the users for 
using too much power and thus generating too much MUI. 

Note that the computation of the optimal power allocations of the users in Algorithm [5] can be performed 
locally and distributively by the users as previously discussed for the Pf NEP, once e^") and 7j are given. 
The computation of 7^ requires an estimate from each user i of the cross-channel between its transmitter 
and the receivers of all SUs being in the coverage radius of user i. This estimate needs to be computed only 
once (before running the algorithm) and updated at the rate of the coherence time of the channel. When the 
computation of 7^ is not possible, one can still use Algorithm [5l setting 7^ = in ()65p . which corresponds 
to solving the optimization problem (j63p with (j){p) = 0, and thus computing just one of the solutions of 
the game ^sisoj Theorem still guarantees convergence of the algorithm, even in the presence of multiple 
equilibria. The lack of any coordination however does not allow to perform any solution selection. 

7.1.1 Efficient computation of the SISO proximal best-response solutions 

In the distributed algorithms described in the previous section, each SU needs to compute the best-response 
solution of his rate maximization problem; see, e.g., (165 p . In this section, we provide an efhcient method for 
computing such a solution. Motivated also by other resource allocation problems, such as [59], we introduce 
next a very general formulation that encompasses the optimization problems studied in this paper, whose 
optimal solution is proved to have a multi-level waterfilling-like expression, and provide an efficient algorithm 
to compute the optimal water-levels (dual variables). 

Consider a parallel additive colored Gaussian-noise channel composed of N subchannels with channel gains 
{^k}k=ii subject to some spectral mask constraints p™'^^ = (p™^^)^-^ and to several weighted average power 
constraints across the subchannels ^k=i^kPk ^ ct, where the inequality has to be intended component-wise 
and Q > 0. Under this setup, the mathematical formulation of the optimization problem in the presence of 
given exogenous prices A = (A^)^-,^ and the proximal regularization in the objective is the following: given 
r > and c = (cfc)^^ 

maximize Y.k=i [log + HkVk) - >^kPk] - § ||p - c||^ 

Ek=i ^kPk < « (68) 
< p < p™^^ 

where we make the following trivial assumptions without loss of generality: < {Hk)^^i < 00; < = 
{wki)^i < 00 for all k = 1,...,N, and linearly independent; < (A^)^-,^ < 00; < {ck)^^i < 00; 
< p"^'^^ < 00; and J2k^i ^kPT"" > 

Problem (|68l) is a convex problem with a polyhedral feasible set; the KKT are then necessary and sufficient 
conditions for the optimality. Solving the KKT system leads to the following structure for the optimal solution 
p*: denoting by fi = > the multipliers associated to the weighted average power constraints, for 
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each k = 1, . . . ,N, we have 



Pk 



Pt 
Pi 
0, 



if \k + ^^' ^k< T+w^ 



pI: Xk + fJ.'^v/k + t{pI- Ck) = i+^Ip* , if J+^^m^sc -rip^^"" - Ck) < Xk + fJ-^^k < Hk + TCk 



if Hk + TCk < Xk + n'^v/k, 



where the multiphers fi have to be chosen to satisfy the complementarity conditions 

N 

< Pi ± ai -^WkiPl>0, i = l,...,Nc. 



(69) 



(70) 



k=l 



The optimal power allocation in (|69p is defined only implicitly, as a solution of a quadratic equation. To 
obtain an explicit and simpler expression, we study first the properties of such equation. 



Lemma 35 Consider the quadratic equation jlk + T {Pk — Ck) = ij^^p^ in the variable pk- The equation has 
the following properties: i) both roots are real, one root is always negative, and one root is always nonnegative; 
a) both roots are decreasing in fik, and Hi) the nonnegative root is given by 



Pk 



Ck 



Hr. 



1 

27 



l^k 



l^k-'T\Ck + 



Hr. 



+ 4r 



(71) 



Proof. See Appendix [HI ■ 
Applying now Lemma [35] to ([69|) [using, in particular, the fact that the positive root of the quadratic 
equation in (I69p is decreasing in the term flk — Xk + /^"^w^], it is not difficult to see that the optimal power 
allocation ([69|) can be equivalently written as 



Pk 



l\Ck 



1 
Ifk 



1 

27 



l^k 



+ 4t 



/c = 1, 



,iV 



(72) 



where each flk — Xk + /^"^w^ and the water-level vector has to be chosen to satisfy the complementarity 
conditions in (I70p . and [x]^ denotes the Euclidean projection of x onto [a, b], i.e., [x]^ — max(a, min(x, 5)). 

It is instructive to verify that the expression in ([72|) particularizes to the multi-level waterfilling when 
r = 0: 



lim 



l\Ck 



1 

Hk 



1 

27 



Afc - 7" ( Cfc + 



+ At 



1 1 

/ifc Hk 



Also, when r goes to infinity, the optimal power allocation is equal to the center of the regularization: 



lim 

T— >00 



Ck 



1 

Hk 



1 

27 



Pk-r\Ck + 



Hk 



+ At 



[Ckfo" 



The numerical computation of the water-level /x in (I72p so that the complementarity conditions in (|70j) 
are satisfied can be done efficiently in practice with a multiple nested bisection method as described in detail 
in Algorithm [6] below. 
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Algorithm 6: Multiple nested bisection algorithm for the computation of the proximal 

best-response in (|72|) . 

(S.O) : Choose some accuracy e. 

(S. 1) : Set = and JI^ = max^ {{Hk - TCk - Xk) /wki}- 

(5.2) : Set = (^^ +7Ii) /2. 

(5.3) : Solve for /i2, • • • 

(S-2 . 1) : Set [£^ = and 7I2 = maxfc | (^Hk - TCk - Afc - Y.j<2 '^kjfJ-j^ /i'^fc2}- 

(S-2 . 2) : Set = {ji^ + ^^2) /2- 
(S-2. 3) : Solve for . . . 

(S-3.1) : Set M3 = and 

7I3 = maxfc I (^Hk - TCk - Afc - Y.j<^ WkjUj^ /"Wfes}- 

(S-3. 2) : Set /X3 = [[£.^+'^3) /2. 
(S-3 . 3) : Solve for fi^, . . . 

(S-3.4) : Using (f72|) . if J2k=i ^ksPt < 0:3 then set 7I3 = fi^, otherwise set fi^ = fi^. 

(S-3 . 5) : If 7Z3 - ^3 > e, then go to (S-3 . 2) . 
(S-2.4) : Using ([72]) . if J2k=i ^fc2Pfc < 02 then set J12 = IJ.2, otherwise set fj.^ = jj,2. 
(S-2 . 5) : If 7Z2 - ^2 > ^' ^^^"^ SO to (S-2 . 2) . 

(5.4) : Using ([72]) . if Ylk=i "^kiPk < o;i then set Jli = fii, otherwise set = fJ-i- 

(5.5) : If - ^ > e then go to (S.2); otherwise STOP. 



The basic idea of the algorithm is to employ a bisection algorithm in /X]^; then, for a given /i^, use a 
bisection algorithm in fi2] and then in /I3 and so on. For the ith bisection level, the interval can be chosen as 
0, maxfc I (^Hk — rck — Xk — X^j<i f^j'^kj^ /wki^ ■ The convergence of the nested bisection method follows 
:'rom the following nontrivial monotonicity result. 

Lemma 36 For a given value o//X]^.,j = [f^i, . . . , fj.^], letp^ il^i-i) denote the optimal power allocation satisfying 
li7^) and the complementarity conditions in \7(j\) associated to the dual variables (Mi-j)- Then, the 

weighted sum power X^^Li WkiP^ monotonic decreasing in fi^ for fixed 

Proof. See Appendix [HI ■ 
We can now state our main result for the computation of the multi- level waterfilling solution in (f72|) . 



Proposition 37 Algorithmic converges in no more than 

n 

i 

iterations, where e is the desired accuracy in the computation of the parameter /i. 



log2 ( max{(i?fc - TCk - Xk) /wki} /e 



Proof. Note that it is implicitly assumed that e is sufficiently small such that the evaluation of the log 
is nonnegative. The convergence of the algorithm follows from the following fact: For fixed the ith. 
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bisection loop converges in at most log2 ^max^ | (^Hk — TCk — Xk — ^j<i t^j^kj^ /e iterations due 

to the monotonicity of X^^^ WkiP'l in Mj, as given in Lemma [36l The worst case number of iterations 

for the ith bisection loop can be easily obtained by taking =0, for all j < i. ■ 

7.2 The MIMO case 

Let us consider now the MIMO games. Because of space limitation, we focus only on the complex NEP ^mimo 
in (jlip : the analysis of the complex NEP ^femto in (|13p follows the same lines, and thus is omitted. 

Given ^mimoj let us assume w.l.o.g. that all the matrices Qj in the game have the same dimensions ut x tit- 
It follows from Proposition [21] (see also Example I25p that Gmimo is equivalent to the complex VI (T^""^"", F'^) 
where F^{Q) = (Ff (Q))^^^ : ^ ^QurxnT^ ^ith each 

Ff(Q) 4 -(VQ,i2,,(Q))* = -Hf ||r„, + J^H,,Q,H,f j H,,. (73) 

According to Proposition l28l the monotonicity/P properties of F'^(Q) on 'P"'^"'° are related to the properties of 
the augmented Jacobian matrix JF'^(Q). To obtain sufficient conditions easy to be checked, let us introduce 
the following I x I matrix T™" G M^^^ obtained by properly "condensing" JF'^(Q) (see Appendix [I]) , and 
defined as (we implicitly assume that all the channel matrices Hjj are full column rank): 



" •vmimo ] 
. -■- pC \ij 



-p(HlfH,fH,,Hn.lNNR,,- if z / j 



if i = J 

(74) 



where denotes the Moore-Penrose pseudoinverse of A, p(A) is the spectral radius of A, and INNRjj is 
defined as 

INNR., ^ ^ , . ^ (75) 

'^least y^rij 

where, in this case, Ajgg^g|- (RnJ is the minimum eigenvalue of the (positive definite) matrix Rn^- Note that 
since Hjj are full column rank, we have h|- = (H^Hjj) ^ ■ 

Using the above definitions. Proposition [24] and Proposition [28l we obtain the following characterization 
for the complex game ^mimo- 

Proposition 38 Given the NEP g^^^o = (-Ri)i=i), the following hold. 

(a) ^mimo is equivalent to the complex VI('P'"^"'°, F*^), which has a nonempty and compact solution set; 

(b) Suppose that matrix Tp^° is positive semidefinite. Then F*^ is monotone on 'P'^^^° ^ therefore, ^mimo is 

a monotone complex NEP; 

(c) Suppose that Tp^° is a P-matrix (positive definite matrix). Then F^ is a uniformly P-function (strongly 

monotone function) on 'p^^^° ■ therefore, ^mimo is a complex P-j- NEP and has a unique NE. 

Proof. See Appendix [I] ■ 
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Corollary 39 The matrix Tp^c"" is a P (or positive definite) matrix if one (or both) of the following two sets 
of conditions are satisfied: for some w = > 0, 



Low received MUI: — ^ wj |p (H|f H^HyHj.) • INNR^A < 1, 
Low generated MUL -1 {p {^f^fj^ij^l) ' INNRy } < 1, Vj 



Vi = 1, 



(76) 



Wj 

It is interesting to observe that conditions for F*' to be a uniformly P-function on 'p^^'^° are the natural 
generalization of those obtained for ^siso to be a P-f game; they thus have the same physical interpretation, 
for which we refer the reader to Sec. 17.11 Based on that, in Proposition [38] we used the same terminology as 
in Definition [TTl namely: ^mimo is a complex P-j- NEP if T™" is a P matrix, whereas is a complex monotone 
NEP if 

Tpc°° is a semidefinite matrix. For these two classes of NEPs we can devise distributed algorithms 
having the same convergence properties and features of those developed in Sec. 15.11 and Sec. 15.21 for real P-j- 
and monotone NEPs, respectively. The main results are briefly listed next; proofs are based on the same 
techniques used to prove Lemma [2l] and Proposition [38l and thus are omitted. 

The case of P-f NEP ^mimo (low-interference regime). When the matrix Tp^° is a P matrix, the unique 
NE of the game can be computed using Algorithm [1] on ^mimoi as stated in the next theorem. 

Theorem 40 Suppose that Gmimo is a complex Px NEP. Then, any sequence {Q(")}^^o generated by Al- 
gorithm m applied to Qmimo converges to the unique NE of the NEP, for any given updating schedule of the 
players satisfying assumptions Al)-A3). 

The algorithm has the same desired features as the one obtained for the SISO case; see Sec. 17.11 

The case of monotone NEP ^mimo (medium-interference regime). When Tpc""" h 0, t/mimo is a monotone 
complex NEP; in the presence of multiple solutions, we need to choose first a merit function quantifying the 
quality of a NE of t/mimo; similarly to the SISO case, we consider the overall interference generated by all the 
SUs: 

i=l j^i 

where WiS are positive given weights. Building on the equivalence between the game ^mimo and the VI ("P""^™, 
F*^) (Lemma [23]), the NE selection problem based on the merit function (f) becomes: 

minimize MQ) 

Q (78) 
subject to Q G SOLc(P"'^"'°, F^). 

A solution of ([78]) can be computed in a distributed way using Algorithm [1] applied to (j78p : it can 
be shown that this corresponds to solving a sequence of perturbed complex NEPs denoted by ^(„) q = 

■pmimo^ i?j(Qj, Q_j) -|- (•, Ilj) -|- — ||» — Qj||^'j V where each player solves the following convex 

optimization problem: given Q_i ^ 0, Qi G C"^", r > 0, e^") > 0, and Ui =Y,j^iWi^^iB.ji, 

maximize Pj(Qj, Q_j) - (Q^, Ilj) - ^ llQi - QJIp. (79) 
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The formal description of the algorithm along with its convergence properties are given in Algorithm [7] and 
Theorem 1411 below, respectively. 



Algorithm 7: NE selection for the complex NEP Qj^ 



(Data) : {e^} i and r > 0. 



(S.O) 
(S.l) 
(S.2) 



Choose any Q^''^ € -pmimo ^^^-^^ ^ center of the regularization Q ^ 0, set e = e^'^^ and n = 0. 
If Q*^") satisfies a suitable termination criterion, STOP. 
For each i = 1, . . . , I, compute Q -"^^^ as 



Q(n+1) 



' Q^eavgma^\Ri(QiM-r'^'"^^)-HQ^,'n^)-^\\Qi-Q^\\l}, UneTi 
Q , otherwise 



(80) 

(5.3) : If Q("+^) is a NE of Ga,e(^\Q^ ^^^^^ update e and the center Q: 

e = e("+i) and Qi = qS"+^^ Vi = l,...,z; (81) 

(5 . 4) : n n + 1 and return to (S.l). 

Theorem 41 Given the optimization problem (|78p . suppose that: i) ^mimo is a complex monotone NEP; ii) 
{e^"'^} is such that e^") — )• and X^^g ^^"^ ~ ^''^^ ''" ^■^ chosen so that T™" + rl is a P matrix. 
Then, the sequence {Q^"^}^o generated by Algorithmic has a limit point and every such point is a solution 

of m- 



A sufficient condition for matrix T™" + rl in Theorem 1411 to be P is 



^ > l^f^j { Y.P (Hlf Hf H,,Hy • INNR,, ) - 1. (82) 



Algorithm [7] has the same desirable properties of Algorithm [5] introduced in Sec. 17. H which we refer to 
for a detailed discussion on its practical implementation. Here, we only observe that, in the scenarios where 
the SUs cannot compute the optimal pricing matrix Hj we can still use Algorithm [7] setting Hj = in (I80p , 
and converge (under assumptions of Theorem I4ip to any one of the solutions of the game ^mimo- 



8 Numerical Results 

In this section, we compare some of the proposed algorithms solving NEP t/siso in © and (/mimo in (111|1 in 
terms of achievable rates and convergence speed. We also compare the performance of our distributed schemes 
with those achievable computing a stationary solution of the related sum-rate optimization problem (under 
the same power and interference constraints). Among all the schemes available in the literature, we consider 
here the pricing-based algorithm proposed in |59| for the SISO case, and the solution methods proposed in 
|591 [60] for the MIMO setting; we slight modified the schemes in |591 160| in order to include the interference 
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constraints we have proposed in this paper. Algorithms in \59\ I60| are the benchmark methods for this kind 
of problems. It is important to remark that algorithms in \59\ [60] are either centralized or require heavy 
signaling and coordination among the users to be implemented. The proposed comparison then sheds light 
on the trade-off between network signaling and system performance (measured in terms of the achievable 
sum-rate). Finally, we contrast our best-response schemes with gradient-response ones |34) . 

Example ^^1 (NE selection vs. no selection). In Figure|H we plot the SUs' sum-rate X]i=i ^i(p) versus 
the iteration index, achieved by the following algorithms applied to the NEP ^siso^ i) The Jacobi version of 
the proximal decomposition algorithm described in Algorithm [2] (green line curve) ; ii) the Jacobi version of 
Algorithm [5] (blue line curve), where i5^(p) is given by (f62|) : iii) the same Algorithm [5] applied to (j63p . where 
(j){p) in (j62p is replaced by —(j)ip) (red line curve); and iv) the Jacobi Dynamic Pricing-based Algorithm 
(DPA) reaching a stationary solution of the sum-rate maximization problem (black line curve) |59) . The 
choice of the merit function —(f)(p) leads to the selection of the NE solution that maximizes the overall MUI 
in the system, which provides a benchmark of the sum-rate variability and an estimate of the worst-case 
performance over the set of the Nash equilibria of the game. Any best-response solution involved in the 
optimization problems is computed using the waterfilling-like expression introduced in Sec. 17.1.11 

We examined the behavior of the above algorithms under the following setup. We considered a CR 
network modeled as a Gaussian parallel interference channel, composed of / = 25 active users and two PUs; 
all the users are randomly placed within an hexagonal cell; the channels of all the links are simulated as FIR 
filter of order L = 10, where each tap is a zero mean complex Gaussian random variable with variance equal 
to 1/L; the channel transfer functions are the EFT of the corresponding impulse responses over = 128 
points (carriers). We focused on two scenarios, namely low interference regime (corresponding to Tg being 
a P matrix) and high interference regime (corresponding to JGiow being positive semidefinite) . The thermal 
noise variance ai^{k) is set to one for all k and i, and the transmit power budget Pi is chosen so that the 
Signal-to-Noise-Ratio (SNR) of each user SNR^ = 10 log^o {Pi/ai'^{k)) = 5dB for ah i and k. The PUs 
impose interference-temperature limit constraints; for the sake of simplicity, we assume the same interference 
thresholds ctj = a 1 for all the SUs, with a = 10"'^ (this choice of a is such that the power budget constraints 
of the SUs are not active at any optimal solution). All the algorithms are initialized by the same starting 
point, chosen randomly in the set J>^^^°^ and are terminated when the Euclidean norm of the error in two 
consecutive iterations becomes smaller than 10"^. In the PTDA, we chose the center p of the regularization 
randomly in "P, r = 3.5, and e^") = e^^'^ /{I + 10 n), where e*^^^ = 0.5; in all the algorithms, the termination 
criterion of the inner loop, if any, is the same as the outer loop. The above choice of the free parameters is 
the result of some preliminary tests; it is important to remark however that the proposed algorithm has been 
observed to be robust against the variation of the aforementioned parameters. 

The following comments are in order from Figured In the case of multiple Nash equilibria, the sum-rate 
performance of the network can vary significantly over the set of the Nash equilibria; the relative sum-rate 
gap between the "worst" and "best" NE is more than 90%. As expected. Algorithm [5] outperforms Algorithm 
[21 which validates the use of criterion (j62p in choosing the NE. Moreover the sum-rate loss with respect to the 
DPA is very limited, and more than acceptable if one takes into account that, to be implemented, the DPA 
requires a significant signaling among the users at each iteration. There are scenarios where such a signaling 
exchange is not feasible, because the users are heterogeneous systems that are not willing to cooperate with 
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each other; in all these cases the proposed algorithm is a good candidate: it outperforms Algorithm [2] while 
keeping the same distributed nature and guaranteeing convergence. When the NE of the game is unique 
[Yg is a P matrix], as expected, both Algorithm [2] and Algorithm [5] converge to the same sum-rate solution. 
Interestingly, this solution seems to coincide also with the one achieved by the DPA. Finally, note that our 
algorithms converge quite fast. 

Ill 

— NE Selection (min MUl), Alg. 5 \ 

— NE Selection (max MUl), Alg. 5 

— No NE Selection, Alg. 2 

— Dynamic Pricing Algorithm (DPA) [59] 



Interference 
^; unique NE) 

"High/medium" Interference 
(monotone NEP; multiple NE) 

/ 

II 
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Figure 4: Comparison of distributed algorithms solving ^siso: Sum-rate of the SUs versus the iteration index. 

Figure [5] shows the average performance of the algorithms i)-iii) considered in Figure U] We plotted the 
average sum-rate versus the SNR = P, with Pi = P and <Tf{k) = 1 for all i and k, achievable at the NE 
reached by the Algorithm [2] and Algorithm [5] The curves are averages over 5000 random channel realization 
chosen so that the JGiow is positive semidefinite. The rest of the parameters are the same as in Figure [4] 
Figure [5] confirms the superior performance of Algorithm [5] with respect to Algorithm [2] that does not perform 
any equilibrium selection. 

Finally, it is worth observing that projection-response algorithms proposed in [34] and |30[ [33] cannot be 
used to solve the P-f NEP ^siso (unless it is also monotone), even if the game has a unique NE. 

Example #2 (Comparison with gradient-response algorithms for monotone Vis). In Figures [6] 
and [7] we compare some algorithms solving the game ^siso under the monotonicity assumption (^siso is a 
monotone NEP); the setup is the same of Figure |4] More specifically, in Figure (H] we plot the SUs rates 
versus the iteration index achieved by Algorithm 2 and the Iterative Tikhonov Algorithm recently proposed 
in |34| for solving monotone Vis. In the latter algorithm we chose the variable step-size sequences 7^ = n~^''^ 
and 6n = n~^''^^ so that (sufficient) conditions for the convergence given in |33| Prop. 15.1] are satisfied (we 
use the same notation as in |34| : see therein for the details). The figure clearly shows that our best-response- 
based scheme converges in a very few iterations, whereas the gradient-response algorithm needs much more 
iterations (two orders of magnitude more) to reach comparable performance. The same convergence properties 
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Figure 5: Average sum-rate of the SUs versus the SNR at the NE of tJsiso- 

as in Figure [6] has been experienced for all the channel realizations we simulated. 

Figure [7] shows an example where Algorithm 1 does not converge because of the high interference among 
the SUs (sufficient conditions in Theorem [33] are not satisfied), whereas Algorithm 2 still converges (in the 
sense of Theorem I19p in a few iterations. Interestingly, we experienced convergence of Algorithm 2 even when 
JGiow is not positive semidefinite, provided that the error sequence {e^")} is properly chosen (e.g., {e^'")} 
goes to zero sufficiently slow), showing that the conditions for the convergence are far from being necessary. 




Figure 6: Typical behavior of gradient-respouse versus best-response algorithms solving the monotone NEP Gsiso- 
rate of three out 25 users versus the iterations achievable by the gradient-response algorithm [M] (red-line curves) and 
the simultaneous best- response algorithm described in Algorithm 1 (blue-line curves). 
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-Best Response Algorithm (Alg. 1) 
-Proximal Decomposition Algorithm (Alg. 2) 
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Figure 7: Comparison of distributed algorithms for the monotone NEP Q^^^°: Rates of the SUs versus 
iterations achievable by GIWFA (Algorithm 1) and PDA (Algorithm 2). 



Example ^3 (NE selection vs. stationary solutions: the MIMO case). We compare here some 
of the proposed algorithms in the MIMO setting. We consider the same scenario as in Figure 1, with the 
only difference that now the transceivers i = are equipped with ht. = nR^ = 3 antennas and there 

are 1 = 5 active SUs. The channels are MIMO frequency-selective (the order of the channels is L = 10 and 
the number of subcarriers is = 128) and are generated in order to guarantee that the matrix Tpc in (l74|l 
is positive semidefinite, resulting thus in a monotone NEP Q'^^'^° , Soft average power shaping interference 
constraints are imposed to each SU along the direction of the primary transmitters (the column of the 
matrices Gpi are the steering vectors representing the directions of the primary receivers p = 1,2 as observed 
from the secondary transmitters i = 1, . . . ,/); all the interference threshold are assumed to be equal and 
set to = 10"^ The best -response of each user cannot be computed in closed form (unless the proximal 
regularization is not included in the objective function), but can be efficiently computed using any nonlinear 
programming solvers (each player's optimization problem is convex). In Figure [8] we plot the sum-rate versus 
the iteration index achieved by the Jacobi version of Algorithm 2 (green-line curve). Algorithm 7 based on 
the merit function ^(Q) defined in ()77p (blue line curve); and iii) the Gauss-Seidel based algorithm proposed 
in |60| to compute stationary solutions of the sum-rate problem (we adapted the algorithm in |60| including 
the interference constraints in the feasible set of the optimization problem) . Figure [8] shows that trade-off 
between performance and signaling that is achievable by the three algorithms: Algorithm 7 implementing 
a NE selection leads to better sum-rates than Algorithm 2 at the cost of almost the same signaling among 
the SUs of classical MIMO IWFAs (a constant price depending on the cross-channel matrices needs to be 
computing by each SU before running the algorithms); higher sum-rates can be achieved using algorithm in 
|60| but at the cost of more signaling among the SUs, which is not admissible in all CR scenarios (note that 
in the MIMO case, the scheme |60) requires the SUs to exchange matrix informations at each iteration of the 
algorithm) . 
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Figure 8: Comparison of distributed algorithms solving the gcLinc ^mimo ■ Suiii-rcLtG of tliG SUs vGrsus iteicitioiis. 

9 Conclusions 

In this paper, we have proposed a novel method based on Vis suitable to study and solve general real 
or complex player-convex NEPs, having (possibly) multiple solutions. The proposed framework has many 
desirable new features, such as: i) it can be applied to real and complex NEPs having no specific structure and 
for which the best-response mapping is not available in closed form or unique; ii) the algorithms proposed for 
computing a NE converge under mild conditions that do not imply the uniqueness of the equilibrium; and iii) 
in the presence of multiple equilibria, one can control the quality of the computed solution by guaranteeing 
convergence to the "best" NE (according to some prescribed criterion), at the cost of limited signaling among 
the players (or no signaling at all). These features make the proposed algorithms applicable to a variety of 
scenarios in different fields; the choice of one scheme with respect to the other will depend on the trade-off 
between signaling and performance that the users are willing to exchange /achieve. The analysis of algorithms 
for complex NEPs hinges on the definition of the VI problem in the complex domain; this new class of Vis 
along with their properties are introduced and studied for the fist time in this paper. 

Finally, to have suitable case studies, we applied the proposed framework to solve some novel NEPs 
modeling various resource allocation problems in SISO/MIMO CR and femtocell systems. The resulting 
distributed best-response algorithms were shown to converge even when current schemes proposed in the 
literature for related problems fail. Numerical results showed the superiority of our approach with respect to 
plain noncooperative solutions as well as good performance with respect to centralized solutions, in favor of 
very limited signaling among the players. 



54 



Appendix 



A State-of-the-art Solution Analysis of Partitioned VI 

We summarize herein the state-of-the-art results deahng with the solution analysis of partitioned Vis, which 
are more general than those given in Theorem [3l by Proposition [71 these results apply readily to player-convex 
NEPs. Note that some of them are drawn from |30| . some others have not appeared in the literature and are 
new in their genre. 

Before stating the main theorem, we introduce the following intermediate definitions. Given the parti- 
tioned VI(Q, F), let F'q'^ denote the normal map of this VI, defined as 

F^-(z) ^ F(nQ(z))+z-nQ(z), z G M", 

where IIq{z) denotes the Euclidean projection of z onto Q. A solution x* of the VI(Q, F) is isolated if there 
exists an open neighborhood J\f of x* such that TV D SOL(Q, F) = {x*}. Lastly, we let deg(Fg"', O) denote 
the degree of the map Fg"' with respect to the bounded open set O, provided that this degree is well defined 
(see |30| Def. 2.1.1] for a formal definition of the degree of a map). 

/ 

Theorem 42 Let Q = Qi with each Qi C M"' being closed and convex {albeit not necessarily bounded). 

1=1 

Let F = (Fi(x))[^-|^ : Q— )-R" be a continuous Pq function on Q. Then, the following statements hold: 

(a) Existence of a solution ([30l Th. 3.5.11]): SOL(Q, F) is nonempty if and only if there exists a 

^rcf ^ {x.f)(^^ G Q such that the set 

L< = I X = (xi)f^i G Q I Fi(x)^(xi - x.f ) < for all isuch that x^ / xf^ } (83) 

is bounded. The "if" part holds true even ifF is not a Pq function |30l Prop. 3.5.1]; 

(b) Boundedness of solutions ([3U1 Th. 5.5.15]).- SOL(Q, F) is nonempty and bounded if and only if 

there exists a bounded open set O containing (Fg''')~^(0) such that degiF^Q^ , O) is nonzero. A sufficient 
condition for this to hold is the existence of a y^*^^ = (xj°^)^^-^ G Q such that the set 

L'< ^ j X = (x,)Li G Q I max F,(x)^( x, - x^^ ) < | (84) 

is bounded. The boundedness of L'^ is a sufficient condition for SOL(Q, F) being nonempty and compact 
even ifF is not a Pq function [301 Prop. 3.5.1]. 

(c) Connectedness of solutions ([30l Th. 3.6.6(a)]): // SOL(Q, F) is nonempty and bounded, then it is 

connected; 

(d) Uniqueness of the solution ([301 Th. 3.6.6(b)]): The VI(Q,F) has a unique solution if and only if it 

has an isolated solution. □ 

Theorem |32] is very general and provides the state-of-the-art conditions for the above properties of a 
solution of a partitioned VI; for instance, Theorem [3] comes as special case. 



55 



We conclude this section providing conditions for the isolatedness of a solution of a the partitioned VI 
(and thus the globally uniqueness); we also present necessary and sufficient conditions for the uniqueness of 
the solution of a monotone nonlinear complementarity problem. The latter conditions are new in their genre. 

Sufficient conditions for the isolatedness of a solution to the partitioned VI can be obtained by introducing 
the critical cone at a solution x* G SOL(Q,F); by definition, this cone is given by: 

C(x^Q, F) ^ J]C(x^Qi,F,) with C(^,Qi,F,) ^ T{Qi; x^) n F,(x*)^, (85) 

i=l 

where T{Qi;x*) denotes the tangent cone of Qj at x* . Exploiting only the Cartesian structure of Q, the 
result below does not require the Pg-property of F. 

Proposition 43 Let F be continuously differentiable on tC. The following two statements hold for a solution 
X* G SOL{Q,¥): 

(a) // 

max (ui )^[ JxFi(x*)u] > for all nonzero u = {u.i)l^i G C(x*;Q, F), 

l<i<I 

then X* is an isolated zero of the VI(Q, F). 

(b) // each Qi is polyhedral and if the complementarity problem: 

C(x*; Q,F) 3 u ± JxFi(x*)u G C(x^ Q, F)* 

has u = as the unique solution, then x* is an isolated zero of the VI(Q,F). □ 

We can obtain sharpened necessary and sufficient conditions for the isolateness, and thus the globally 
uniqueness, of a solution to a monotone NCP(F): < x _L F(x) > 0, where F is monotone on M" in the 
usual sense, i.e., 

(x-y)^(F(x)-F(y)) > 0, Vx,yGM;. 
To state this condition, define the following index sets associated with a given solution x* of this NCP: 

— {j I > = -^j(x*)} (strongly active constraints) 
/5* = {j \ X* = = Fj(x*)} (degenerate constraints) 
7* — {j \ X* = < Fj{'x*)} (inactive constraints). 

The first two parts of the result below provide necessary and sufficient conditions for x* to be an isolated 
solution of the monotone NCP(F). In particular, part (b) of the result shows that by solving finitely many 
convex programs, each of which is the minimization problem of a linear objective function over a convex set, 
one can definitively conclude whether a given solution to a monotone NCP is isolated or not. Unlike parts 
(a) and (b), part (c) provides a matrix-theoretic sufficient condition for the isolatedness of the solution. 

Proposition 44 Let F be a monotone function on M" and let x* be a solution of the NCP(F). The following 
two statements hold: 
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(a) X* is an isolated, thus unique, solution of the NCP if and only if 

F„^(x„^,X/3^,0) = 
<^I3^ ± F^^(X„,,X;3^,0) > 



/3* 



(b) IfF is continuously differentiahle near li* and Jxq^Fq,^(x*) is invertihle, then x* is an isolated solution 

of the NCP(F) if and only if for every index subset C with complement (3 = f3^ \ P, the convex 
program 

minimize > xj. 

subject to Fq,^(xq,^,X/3^,0) = ^.g^^ 

F-g(xQ,^,X/3^,0) = < x^ 

and F^(xq^,x^^,0) > = x^ 

has x-g = as the unique solution. 

(c) // F is continuously differentiable near x* and 

Jx^FgCx'^) Jx.F^Cx-^) 

is invertible for every index subset [3 C (3^, then x* is an isolated solution of the NCP(F). 

Proof. To prove (a), suppose that x* is the unique solution of the NCP(F). Let (xq,^,x^ ) satisfy the 
left-hand side of (I86p . Since the mapping 

Fa, (xq,, , X/3^ , 0) 

F/3,(x„,,X/3^,0) 

remains monotone, it follows that for all r > sufficiently small, (x*^ + r(xQ, — x*^), rx^^, 0) is a solution 
of the NCP(F). Hence the right-hand side of (|86p holds. Conversely, suppose x is an alternate solution of 
the NCP(-F). Without loss of generality, we may assume that x is sufficiently near x* so that Xq, > and 
F^^(x) > 0. Hence F„, (x) = and x^ =0. By ([86|l . we deduce that x = x*. This completes the proof of 
part (a). 

For part (b), we first show that the feasible region of ()87p is indeed convex. Let (x^^,a;^ ) for £ = 1,2 
be two feasible solutions to ([87]); these solutions must satisfy the left-hand conditions of (I86p . In addition to 
(j88l) . the map 

\ / Fq,, (xq, , x-g, 0) 

/ V F^(Xa,,X3,0) 

is also monotone; thus for all r € [0, 1], z(r) = r(x^^,x^ ,0) -|- (1 — r)(x^^,x^ ,0) satisfies F^, (z(r)) = 0, 
F-g(z(r)) = < z-g(r), and Fg(z(r)) > = z^(r). This establishes the convexity of the constraint set of 
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(|87|) . By the implicit function theorem, the nonsingularity Jx„^Fci^(x*) imphes that x*^ is an isolated, thus 
unique, zero of the equation Fa^{xa^, 0, 0) = 0. It is now easy to deduce that part (b) follows from part (a). 
□ 



B Proof of Theorem [H] 

It is sufficient to show that, under the assumptions of the theorem, the best-response mapping is a block- 
contraction (i.e., a contraction with respect to some block- maximum norm); the latter property indeed 
guarantees that conditions of the asynchronous convergence theorem in |531 Prop. 2.1] are satisfied by the 
asynchronous best-response algorithm described in Algorithm [TJ 

We introduce first the following norms: Given the set of nonsingular matrices C^'s coming from (I23p . the 
block-maximum norm on M", with n = ni -\- . . . + nj, is defined as 



I 11^ — 
X Ki^^i, — max 



for X = (xi)f=i G M", (89) 



where ||»|| is the Euclidean norm on R"' and c = [ci, . . . , c/]'^ > is any positive weight vector; the (weighted) 
maximum norm on is defined as (see, e.g., [53) ) 

ll^ll^vec- max for x G R-^; (90) 

' i=l,...,I Ci 

and the matrix norm ||-||^ j^j^^ on M.^^^ induced by |HI^ vec ' given by 

l|AC,^at = max-J]|[A],,|c„ for A G M^x^ (91) 

Under the under the P property of Tf, let us introduce the best-response mapping 

i3(x) ^ (i3i(x_i))f=i : Q 9 X ^ Q, with i3i(x_i) ^ avgmm^^^Q j,{^i, x.^). (92) 

The contraction properties of ^(x) are given in the following, where Ff is defined in (f22|) [note that the P 
property of Tf, implies the strong convexity of each /j(-,x_j) for any x„j G Q_i, and thus a™'° > for all 
i; hence, Ff is well-defined]. 

Proposition 45 Suppose that Tf in \2S^) is a P matrix. Then, the best-response mapping B{x) is a block- 
contraction, i.e., there there exists some c > such that 

ll^(x) - ^(y)||^i„,k < «c ||x - y||Ci„,k Vx,y G Q (93) 

with ae = ||FF||^,inat < 1- 

Proof. For any two vectors x, y G Q we have by the minimum principle 

{zi - Bi{x)f V^Ji (e,(x), x_,) > Vz, G Qi, i = l,..., I, 
{zi - B,{y)f V^Ji (e,(y), y_0 > Vz,, G Q^, i = 1, . . . , /. 
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Substituting Zj = Bi{y) into the former inequality and Zj = Bi{x) into the latter, adding the two resulting 
inequalities we obtain with Zj = tj(;Bj(y), y_j) + (1 — tj)(;Bj(x), x_j) and some ti G (0, 1): 



< (5i(x) - ^^(y))^ (VxJ. (e,(y), y_,) - V^J, (^^(x), x„i)) 

= (fi,(x) - B.(y))^ V^^,^/, (z) (fi,(y) - ^,(x)) + (^,(x) - ^.(y))^5^ V2^,^./, (z) (y, - x,) (95) 

where the equality follows from the application of the main-value theorem to the univariate, differentiable, 
scalar function 

[0, l]Bti^ {Bi{x)-B^{y)fV^Ji{ti{Biiy),y^^) + {l-ti){B^{x),x^,)). (96) 
Using the definition of a™"^ and /S™*^^ as given in (j23p , we deduce from the inequality in ()95p 

||^,(x) - B,iy)\U af^ < ||x, - y, ||^. , i = !,...,!, (97) 

(the inequality in (j95p is trivially satisfied if ||;Bj(x) — Hj(y)|| • = 0). Introducing the vectors eg = {ejs.)^^^ 
and e = {ei)^^^ with eg^ = ||'Bj(x) — Sj(y)||,j and = ||xj — yj||j, using the definition of Tp in (f22|) . and the 
fact that a™™ > for all i, (j97p can be written in vectorial form as 

eB<rFe, Vx,yGQ. (98) 

It follows from (f98|) that, for any given c > 0, we have 

ll^(x) - S(y)|lblock = l|eB||^,vec < l|rF||^,^at l|e||^,vec = l|rF||^,^at 11^ " yllblock (99) 

which proves the inequality in (|93p . To complete the proof we need to show that IITfII^ j^^g^^ < 1 for some 
c > 0. Invoking |5H Lemma 13.14] and \53\ Cor. 6.1], we obtain the desired result: 

Tf is a P-matrix <^ 3c >0 such that ||rF||^^jnat < 1- (100) 



C Proof of Proposition [TT] 

Since ^ is a monotone NEP, the VI associated with the NEP ^^,y-the VI(Q,F + r(I - y)) with F = 
(Vxi/j)f=i~is strongly monotone on Q; it follows by Theorem [8)^d) that Q-r^y has a unique NE for any given 
T > and y G M". Let us denote such a unique NE by S^(y) = SOL(Q, F + t(I - y)). 

Necessity: Let x* G Q be a NE of the monotone NEP Q. By Proposition [71 x* G SOL(Q, F); then the 
following hold: 

(x - x*)'^ F(x*) > Vx G Q 
4^ (x - x*)^ (F(x*) + r(x* - X*)) > Vx G Q 
^ X* = SOL (Q, F + t(I - X*)) = S,(x*), 
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implying that x* is the unique NE of Gt,x*- 

Sufficiency: Let x* be a NE of ^r,x*- Then, we have x* = St-(x*), which leads to the desired result: 



(x - S,(x*))^ (F (S,(x*)) + r(S,(x*) - x*)) > Vx G Q 

^ (x-x*)^F(x*) > VxgQ. 



D Proof of Theorem [22 



To prove the theorem we hinge on the theory of Vis. We preliminary observe that the game G is equivalent 
to the VI(Q,F), with F = (Vxji)i=i (Proposition [7]); SOL(Q,f) is thus also the solution set of the VI, i.e., 
SOL(Q, F) = SOL(Q, f). Moreover, still invoking Proposition [71 we have that the game ^(n) ^(n) in Step 2 
of Algorithm m is equivalent to the VI(Q,F("')), where 

F(")(x) ^ F(x) + e(")V(^(x) + r (x - x^")). (101) 

Observe that, under the assumptions of the theorem, F^") is strongly monotone [Definition [2]^iii)], implying 
that the VI(Q, F^")) has a unique solution [Theorem [3]^c)] and thus x^""''^) in Step 2 of Algorithm [H is well 
defined at each iteration. Moreover, denoting by S the solution set of (j33p . assumptions of the theorem, 
ensure that S is nonempty, bounded, and convex. Let us introduce for each n, 

A idist(x("),5) = ^ ||x(") -P5(xW)||2, 

where Ps{y) — ai'gniinxg^Hx — y|| denotes the Euclidean projection on the nonempty, closed, and convex set 
S. Then, to prove the theorem it suffices to show that the sequence {J*-"^} converges to zero. Observe first 
that, since x'^""''"^) at Step 2 is the solution of the game G^ ^(rt) ^{n)—the VI(/C, F^"))— we get, for any y G Q, 

F(x('^+i))+eWV(/.(x("+i))j^(y-x("+i)) > r ( x^ - x('^+i) )'^( y - x("+i) ). (102) 
We can write 

Sin+l) _ §in) ^ 1 II ^{n+l) _ p^(x(«+l)) ||2 _ ^ || _ P5(xW) f 

< ^ ||x("+l) - P5(x("))f - ^ II xW - P5(x(")) f 
= _i||x("+i) - x(")|p + (x^*") - x("+i))'^(P5(x(")) - x("+i)) 



< - - ||x("+i) - xW||2 + - F(x("+i))^(P5(x(")) - x("+i)) 
2 r 



(103) 



An) 

+ L_ V</.(x("+i))^(P5(x(")) - x("+i)) 



r 



(c) 1 e-(") 

< -- ||x("+i) - x(")||2 + — V(/.(x("+i))^ (P5(x(")) - x("+i)) 
2 r ^ ^ ■ 
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where: (a) follows readily from the definition of projection; (b) comes from (fT02ll evaluated at y = P5(x(")) € 
Q; and (c) can be obtained by observing that since P5(x(")) G 5 C SOL(Q,F) and x("+i) G Q, we have 
F(P5(x(")))^ (x("+^) - P5(xW)) > 0, which yields in turn, by the monotonicity of F [Definition E^i)], 
F(x("+-'^))^ (P5(x(")) — x("+^)) < 0. We now distinguish three cases. 

Case 1 : Eventually, < 0. 

In this case the nonnegative sequence {(5(")} is (eventually) non-increasing and therefore convergent. Let us 
denote by no the index from which all are non positive, and let us consider n > uq without loss of 

generality. Since S is bounded, this implies that also {x*^"-'} is bounded. Furthermore, it follows from (|lQ3p 
that {<5{"+i) - <5(")} converges to zero and (5^"+^) - < -(1/2) ||x("+i) - x^^^f, which shows that 

lim ||x("+i) -x(")|| = 0. (104) 

n-->oo 

Summing (jl03p from no to n — 1, we get 

n— 1 
i=no 

bmce converges and V^^'^ < 0, assumption (ii) in the theorem implies that lim sup = 0. Then, 

there exists a subsequence J such that 

lim = 0. (105) 

nSJ 
n — ^oc 

Since {x^")} is bounded we may assume, without loss of generality, that lim ns./ x^*^) = x. Note that, since 

n— J-oc 

Q is closed, x € Q. We show that actually x € SOL(Q, F). If this is not so, there exists a point y G Q such 
that F(x)^(y — x) < 0. Since x^") is the solution of the VI (/C, F^")) in Step 2 of the algorithm, we can write, 

F(x(")) + t(x(") -x("-i))l^(y-x(")) + e("-^)V0fx("))^(y-x(")) > 0. (106) 



By continuity, the definition of y, the boundedness of {x^")}, and (1104^ . we have, without loss of generality 
(after a suitable renumeration) , 



lim F(x("))^(y-x(")) < 0, (107) 

Tie J 

n—^oo 

lim t(x(") - x("-^))^(y - x(")) = 0, (108) 

jiG .7 

n— foo 

lim e("-i)V0(x("))^(y - x^") ) = 0, (109) 

n—^OG 

which contradicts (fT06]l . Therefore x G SOL(Q,F). 
Thanks to (fTOD we have lim x^""^) : 

n£j,n—^oo 

x) = 0. But the convexity of (p implies that 



Thanks to (jlQ4p we have lim x^" = x. Therefore, by (jlOSp and continuity, we get V(j){x)^{Ps('x)- 

n£j,n—^oo 



(PiPsiSi)) > <^(X) + VcPiScY (PsiSc) - X) = (PiSc 
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thus showing that x E 5. Therefore we get hm nej = 0. But since the whole sequence IS 

n— )-oo 

convergent, this implies that the entire sequence converges to 0, thus concluding the analysis of Case 1. 

Case 2: The two index sets J and J are both infinite, where J = | n | V^^^ ^ } '''^^ 

|n G J I -i||xW-x(-"i)f + £^yW > o|. 

By (fT03D . if n G J it might happen that J^") > (J^""^), while if n J then necessarily 5^") < ^("-i). 
Therefore, since J is infinite, to prove that {(5(")} goes to zero it is enough to show that the subsequence 
{(J*-"^}/ converges to zero. To this end, first observe that for every n G J it holds that 

^{n-l)y(n) > ^||xW-x("-l)||2. (110) 

The sequence {x^^^jj is bounded since the definition of V^'^\ (|11U|) and convexity imply (j){Ps{'x.^'^~'^^)) > 
i^(x(")). But 0(^5 (x^"^-'^))) is the optimal value of (j33p and therefore is a number, say /?, that does not depend 
on the iteration n. Therefore, since x^""^) belongs to Q, we have that {x("^}j is bounded. By continuity, 
also {V^^^}j is bounded. Hence, since {e^"^} converges to 0, (|110p yields 

lim llx^'^) -x("-i)|| = 0. (Ill) 



n — J-oo 



Since {x^'^^jj is bounded, it has limit points. Let J C J be a subsequence such that lim „g j = x. Reasoning 

n — >QO 

exactly as in Case 1, (the only difference is that instead of (|104p we use (jllip ). we may deduce that x G 
SOL(Q,F). By continuity, V(/)(x)^(P5(x) - x) > 0. Thus x G 5; hence lim „gj (^^^^ = 0. Since this 

reasoning can be repeated for every convergent subsequence of {x("^}j, we conclude that lim „gj 6^"^^ = 0, 

n — J-oo 

thus concluding the analysis of this case. 

Case 3: The index set J is infinite while J is finite. In this case, the sequence IS non-mcreasmg 

eventually. Therefore {(^^"^} converges, implying that {x(")} is bounded, - <^(")} converges to zero 

and therefore, by (|103p . also (|104p holds. At this point, we can proceed exactly as in Case 1 and Case 2 to 
prove that converges to zero, thus concluding the proof of the theorem. □ 



E Proof of Lemma [24 



To prove the lemma it is sufficient to show that the first-order Taylor expansion as given in (H2l) holds for 
the function /; the rest of the proof follows similar steps as those used to prove the minimum principle for 
(real- valued) functions of real variables and thus is omitted; see e.g., |58| . 

Before proving the lemma, it is useful to introduce a real-coordinate representation of real- valued functions 
of complex matrices and establish the connection between standard derivatives of this representation and the 
M-derivatives of the original functions of complex variables. 

The complex space C"^™ of dimension n • m has a natural structure of a real space M^nm dimensions 
2n ■ m; this comes readily, e.g., from the following isomorphic transformation: 

vec (Re(Z)) 
vec(Im(Z)) 



- A 



l)2nm 



(112) 
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n 



z G 



n2nm 



^ A 
Z = 



for some 7j ^ Z 



;ii3) 



For the sake of simplicity, in the following, we will denote hy Z = C"^"* the original complex space and by 
Z the elements of ^; 7^ = M^nm ^-j^ ^j^g . m-dimensional space of real vectors in the form z, i.e., 

Z.R ^ vec(Re(Z)) 
z/ vec(Im(Z)) 

elements of IZ will be denoted by z, and partitioned as in (jll3p . 

Given a real- valued function of complex matrices / : 2 — )■ M, the representation of /(Z) under the iso- 
morphic transformation (jll2p is denoted by /(z) = /(Z). Note that if /(Z) is M- (continuously) differentiable 
on Z then /(z) is (continuously) differentiable on TZ. Moreover, we can easily establish the connection be- 
tween the Jacobian of /(z) and the Jacobian and conjugate Jacobian of /(Z) [cf. (j39p ]. as shown next. By 
definition, for any z £ TZ, the Jacobian of /(z) is 

-df{z) dfiz) 



dzl 



dzj 



Dj{z) ^ {VJiz)f 
Using ([37|) and ()39p . it is not difficult to see that, for any given Z B Z 

DzfiZ) 



[DUi^), D,J{z)] . 

> z £ TZ, the following hold 



^z*/(Z) 



9/(Z) 



1 



[D^Jiz)+j-D^J{z) 



(114) 

(115) 
(116) 



f?vec (Z* 

which provides the desired relationship between Dzf{z) and Dzf (Z) and Dz*f (Z). 

Exploring the above equivalences, we can now readily prove Lemma [23] leveraging on standard real calculus 
results. Given a real- valued convex and continuously M-differentiable function / : /C — )• M on /C, the first-order 
Taylor expansion of /(Z) = /(z) at /C 9 Zo(<J^ zq) exists and it is given by: 



/(Zo + AZ) - /(Zo) = /(zo + Az) - /(zo) 
~ DJizo) ■ Az 

= Re { [D,J{zo) - j ■ D,J{zo)] [(Az)^ + j • (Az),]} 
= 2Re{L»z/(Zo)vec(AZ)} 

= 2Re{tr((Vz/(Zo))^Az)} 
= 2Re{tr((Vz./(Zo))^Az)} 
= 2 (AZ, Vz./(Zo)) 



(117) 
(118) 
(119) 
(120) 
(121) 
(122) 



where ()117p follows from the first-order Taylor expansion of real-valued functions of real vectors (see, e.g., |58|): 
pl9ll follows from (fTTSl) : (fT20]l is due to Dzf (Zq) = vec (Vz/(Zo))'^ and the property vec (A)'^ vec (B) = 
tr (A'^B) for any A, B e C"^*"; and in (fT2T]l we used the fact that / is real and thus (Vz/(Zo))* = Vz*/(Zo). 
This completes the proof. □ 



F Complex Matrix Derivatives in Example [25 

We derive here the expressions of the (conjugate) derivatives used in the Example [25] and Example 125 [revisited . 
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In order to obtain the expression of the augmented Hessian "HzZ'/lZ), we need to compute V|2*/(Z) 
= Dz (Vz*/(Z)) and V|.z*/(Z) = Dz- (Vz-/(Z)). We prehminary compute Vz/(Z) and Vz*/(Z). 

RecaUing that gOl Prop. 3.12] dlndet (Z) = Tr (Z^MZ) for all Z such that detZ / 0, with dlndet (Z) 
being the (complex) differential of Indet (Z), we have (up to a constant positive factor) 

df{Z) = Tr(^(R„ + HZH^)"^HdZH^j =vec'^^(^H^(R„ + HZH^)"^H)^^ vec((iZ) (123) 

which, using the identification rule as given in [40^ Table 3.3], leads to the following Jacobian matrices of /: 



L>z/(Z) = vec^ ^(^H^ (R„ + HZH^) ^h)^^ and L>z-/(Z) = 0, 



and thus [cf. (l39 



Vz/(Z) = (h^ (R,, + HZH^) ^ h) and Vz*/(Z) = 0. 
It follows from (fT25]l that 

Vz*/(Z) = Vz./(Z) + Vz./(Z)* = (Vz/(Z))* = H^ (R„ + HZ^H^)"' H. 

Given (I126p . we can now compute V|;2*/(Z) and V|,2*/(Z)- The differential of Vz*/(Z) is: 



vec 



d Vz./(Z*) 



vec H^d(R„ + HZ^H^) H 

- vec [G(Z)(dZ)-^G(Z)] 

^ V ' 

G(Z)^H«"(R„+HZ»H«)"^H 

- [G{Zf G(Z)] vec [{dZ*f] 

- [G{Zf G(Z)] K„2„2vec [{dZ* 



(124) 



(125) 



(126) 



(127) 
(128) 

(129) 
(130) 



where in (fT28]l we used the rule dZ~^ = -Z'^{dZ)Z~^ [40l Prop. 3.8]; (fT29]l follows from the property 
vec(ABC) = (C^ ® A) vec(B) |40| Lemma 2.11]; and in the last equality we introduced the commutation 
matrix K„2„2, which is the v? x r? permutation matrix such that vec(A"^) = K„2„2vec(A) |40| Def. 1.8]. 
It follows from (|13Up and the identification rule |4m Table 3.3] that 



V|.z./(Z) 
V|zJ(Z) 



[G(Z)^ ® G(Z)] K„2„2 



(131) 
(132) 



which leads to the expression of the augmented Hessian 'Hzz*/(Z) as given in 



G Proof of Propositions [28] and [29] 

It is sufficient to prove only Proposition 1291 Proposition [28] is just a special case. To do that, we need the 
following intermediate result. 
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G.l Mean- value theorem for functions of complex variables 



We provide here a version of the mean-value theorem that is suitable for real-valued functions of complex 
matrices. We focus directly on the specific function that we need to prove Proposition [29l 

Given a continuously M-differentiable matrix function F*^ : K, — t- C"^'" on the convex and closed set 
fC C (C">^™ and a point AY € (j^"x™^ us consider the real- valued function of complex matrix variables 



For every two points Zi, Z2 G /C, with AZ 
defined as: with Z{t) = Zi + tAZ, 



For some t G (0, 1), we have 



5(Z) 4 (AY, fC(Z^ 
= Z2 - Zi, let h{t) : [0,1] - 

[0, i]3t^ h{t) ^ (ay, F<C (Z(t)) 



(133) 

be the real-valued scalar function, 

(134) 



g{Z2) - <7(Zi) = h{l) - h{0) = h'{t) (135) 

where h' (t) is the first order derivative of h(t) [note that h is continuously differentiable on (0, 1)], and the 
last equality in (jl35p follows from the mean- value theorem applied to the function h{t). To compute h' (t) we 
use the chain rule for complex matrix derivatives |40) as shown next. Rewriting h{t) as 



- tr f AY^ F"^ {Z{t))) + - tr f AY^ F^ (Z(t))* 



and using 



h{t) 



Dz tr (^AY^ F^ [Z] 
Dz' tr (^AY" F^ (Z; 

i:>ztr (ay^f'^(z)' 

Dz* tr (aY^ F^(Z)* 



(136) 



vec(AY*)^ I?zF^(Z) 
vec(AY*)^L>z*F^ (Z) 
Dz* tr (^AY^ F^ (Z 
Dz tr f AY^ F*^ (Z) 



we have 



h (t) = Dzit)h{t) DtZ{t) + Dz^tyHt) DtZ{ty 

= ^ vec (AY*)^ DzF^ (Z(t)) vec (AZ) + ^ vec (AY)^ DzF^ {Z{t))* vec (AZ) 



(137) 



+ ^ vec (AY*)^ Dz*F^ {Z{t)) vec (AZ*) + ^ vec (AY)^ Dz*F^ (Z(t))* vec (AZ*) 



^ vec ([AY, AY*])^ 



vec([AZ, AZ*]), 



(138) 



Dz¥^{Z{t)) Z)z*FC(Z(t)) 
DzF^{Z{t)Y Z)z*F^^(Z(t)) 

where in ()137p we used the chain rule |40| Th. 3.1]. Using (|138p and the augmented Jacobian matrix JF^(Z 
introduced in (j46p . we can rewrite (|135p in a compact form as 



<7(Z2)-5(Zi) = vec ([AY, AY*])^ JF^(Z(t)) vec ([AZ, AZ*]). 
which is the desired result. 



(139) 
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G.2 Proof of Proposition | 

We prove only (a)-(c); the proof of (d)-(e) follows similar steps. 

Sufficiency part. For (a)-(c), it is enough to prove only (c). Given two points Zi,Z2 G /C, let us define 
AZ = Z2 — Zi; we have, for some i G (0, 1), 

'Z2-Z1, F'C(Z2)-F'^(Zi)\ =vec([AZ, AZ*])^ JF'^(Z(t)) vec ([AZ, AZ*]) , (140) 



where the equality follows from (|139p . Since Zi, Z2 G /C, we have that AZ G S/c] moreover Z(t) G /C (due to 
the convexity of /C). It follows from (I140p that if there exists a constant Csm such that vec ([Y, Y*])''^ JF'^(Z) 
vec([Y, Y*]) > {csm/2) \\Y\\p for ah Y G cS/c and Z G /C, then F^ (Z) is strongly monotone on /C. 
Necessity part. Let us focus on the strongly monotonicity property only; monotonicity is obtained in a 
similar way. Suppose that F*^ (Z) is strongly monotone on /C with constant Cgm > 0. Let us show first that 
vec([Y, Y*])^ JF'^(Z)vec([Y,Y*]) > (csm/2) ||Y||^ for aU Y G 5^; and Z G ri(/C), where ri(/C) denotes the 
relative interior of /C (see |58| Ch. 1.4] for the definition of ri(/C) and its main properties). Then, we have 



vec([Y,Y*])^JF^(Z)vec([Y,Y*]) = - vec ([Y, Y*])-^ lim - 

2 t4-0 t 



vec (F^(Z + tY) -F'C(Z)) 
vec (FC(Z + tY) -F'C(Z))* 



(141) 



= ^ lim 1 (tY, fC(Z + tY) - F^(Z)) (142) 

> ^ lim^ ||tY||^ = ^ ||Y||^, VYG5ycandZGri(/C), (143) 

where the equality in (|14ip follows from the (M-) differentiability of F*^ (Z) on /C [ ()14ip can be proved using 
the same approach as in the proof of Lemma [2H but applied to vec (F*^(Z))]; in ()142p we used the definition 
of inner product (|43p ; and in (jl43p we used i) the fact that for sufficiently small i > 0, Z + tY G /C [since 
Z G ri(/C)], and ii) the strongly monotonicity of F'^(Z) on /C. 

Next, let Z G /C but Z ^ ri(/C); by |58| Proposition 1.4.1(a)] there exists a sequence {Z^} C ri(/C) such 
that Zfc —7- Z. By (|143p evaluated in each Z^ G ri(/C) and the continuity of JF'^(Z), it is follows that 

vec([Y,Y*])^JFC(Z)vec([Y,Y*]) = lim vec ([Y, Y*])^ JF^(Zfc)vec ([Y, Y*]) > ^ ||Y||^ , 

for all Y G Sic- This completes the proof of the necessity part. 

Note that Proposition [29l reduces to Proposition [28] if the set /C has nonempy interior. Indeed, when this 
happens, Aff(/C) = C'^^''" and thus 5/c = C'^^'". 



H Proof of Lemma [35] and Lemma [36 
H.l Proof of Lemma 



Consider the quadratic equation /i^ + r {pk — c^) = ij^^p^ or, equivalently, 

rHkpi + {t + Hk ifj-k - TCk)) Pk + (Afc - TCk - Hk) = 0. 
We can write the equation as ap\ + hpk + c = with the obvious identifications. 



66 



i) To show that both roots ^—6 it — Aacj /2a are real, it suffices to show that 6^ > 4ac. After several 

manipulations, it follows that 

62 - 4ac = [r - (h " rck)f + AtHI > 0. 

ii) We will show the decreasingness through the derivative: 



(up to a nonnegative scaling factor) = —Hk\/b'^ — 4ac ± if^ (r — -ff^ (/U^ — rcfc)) 



-//fc V [r - //fc (Afc - TCk)r + 4ri?^ ± (r - //fc (/ifc - rcfc)) 



< -^fc (^V ["^ - -H'fc (Afc - ^Ck)f ±{T-Hk iflk - TCk)) 

= -Hk (|r - Hk ifj-k - TCk)\ ±{t- Hk (/Ufc - rcfc))) 

< 0. 



iii) It suffices to show that b + yjb'^ — 4ac is positive for /i^ = combined with the decreasingness. Indeed, 
for jlj^ = 0, b = T — HkTCk and 

b'^-iac = [T + HkTCkf + 4tHI 

= [t - HkTCkf + 4ri7fc (rcfc + Hk) . 

It follows that 



b+Vb"^- Aac = {t - HkTCk) + \I[t - HkTCkf + At Hk {TCk + Hk) 

> {t - HkTCk) + \t - HkTCk\ 

> 

where we have used c^; > 0. 

iv) The expression for the positive root follows from i^—b + y/b"^ — Aac^ /2a together with 6^ — Aac = 
[t - Hk (/ifc - TCk)f + AtHI, b = T + Hk{fj,k- TCk)-, and a = t Hk- 

H.2 Proof of Lemma 1361 

To simplify the exposition, we will absorb the fixed parameters /^i-j^i in \k as <— Xk + fJ'i-i^i^k,i:i-i7 we 
will denote and Wki simply as fi and Vk, and we will denote /Xj_|_i:jv^ and Wk^i+i-.N^ as p and w^. Therefore, 
instead of writing + fJ.'^'Wk we will more explicitly and conveniently write A^ + fiVk + p^Wk- With this 
new notation, the goal of the proof is to show the decreasingness of '^k=i '^kPk (a*) — 0- 

First of all, observe that, for fixed parameters pi.^_i, p* {jj) is unique and continuous on [0,oo), as shown 
next. Given ^, the vector p* (^u) defined in (j7ip can be interpreted as a solution of the VI(/C, Vp£(-, p*)+/uv), 
where /C^{pGR^:0<p< p--} and C{-,p) ^ - Zk=i N (1 + ^kPk) - {Xk + P^y^k)Pk] +i ||p - cf , 
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with p* satisfying the complementarity conditions < p* _L '}2!k=i^ kPkil^) ^ 0- It is not difficult to see 
that Vp£(-,p) (and thus Vp£(-,p*) + /iv) is strongly monotone on /C, for any given p. It follows that: i) 
the solution p* (/x) of the VI(/C, Vp£(-, p*) + /iv) is unique, and ii) p* (— /i) is a co-coercive function of jj, [30^ 
Prop. 2.3.11] implying that p* (/x) is Lipschitz continuous. 

Now, the optimal solution p* (//) in (|72p contains three cases depending on whether = 0, = p™*^, 
or p^ € (0, p^'^^). The parameter /x naturally partitions [0,oo) into intervals where one of the three cases 
is active for each k. It suffices to assume that ^ lies in the interior of one of these intervals and show then 
the decreasingness of X^^x ^kPk (a*) the interior of each interval. This, together with the continuity of 
'Ylik=i '^kPk (a*)' implies the decreasingness over the whole interval [0, oo). 

Depending on the value of /x, some of the inequalities in Xl^i ^kPX — ^ ^^^^ active and some inactive. 
The inactive ones will have the corresponding elements in p equal to zero. From now on, we will only consider 
the active ones and, for simplicity of notation, we will assume that only the active ones are included in the 
vectors and a so that we can effectively write X^^^ '^fcPfc — (Note that if no constraint is active, then 
"^^=1 "^kPk (a*) trivially decreasing in /i since each is decreasing in ji.) 

Since, by assumption, /x is in the interior of an interval where the three cases of the power expression in 
()72p do not change for each fc, it follows that in a sufficiently small neighborhood of /i the indexes k where 
either p^. = or p'^ = p^'"^^ do not change and, therefore, do not affect any of the terms in X^^^ ^kPk — ^■ 
We can then ignore those indexes and only consider the ones where p^ G (0, p™^^), i.e.. 



f^k 



""'^Ifk 



+ At 



where fii. = + /UUfc + P^^k- At this point, we can write the system of equations as ^kPk ~ where 
the summation is over the indexes for which p^ € (0, p™**^). There is yet another simplification we can 
make: if the vectors are not linearly independent, we can just use another set of vectors with reduced 
dimension represent the same system of linear equations. With that reduction in mind, we will still use the 
same notation and write Xlfc^fcPfc = ^ 

The active constraints can then be rewritten as the following system of equations: 



f (a^, p) 




fJ-k 



l^k 



+ 4r 



Q = 



where the summation is only over the indexes that satisfy G (0, p^^^) and a contains a minus the effect of 
the indexes that we are not considering in the summation (i.e., those for which either p^. = or p^ = p^^^). 

At this point, we can invoke the implicit function theorem to characterize p (^u). In particular, the partial 
derivatives of f (/i, p) are 



5f 

_di_ 



4^ ^djlkdp^ 



dpi 



^kVk 
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and then 



E 



-1 



where the inverse exists because of the assumption of w/^. being Unearly independent and dp^/dfif. < for all 
k such that G (0, p^^^) (Lemma [35]). 

We can now write the derivative of the weighted sum power as 




dp jp) 



^ ^Pk 



T 



si!"'""- 



Proving that this derivative is nonpositive is equivalent to showing (since dp'^/dfij^ < 0) 

\ -1 / 



dpi 



Vk^k 




dpi 



E 



dpt 



^kVk 



< 



dpi 



dh 



(144) 



which follows from Cauchy-Schwartz's inequality for random vectors: let x G M" and y G be two random 



vectors such that E 



X 



< oo, 



< oo and E [yy"'"] is nonsingular, then E [xy-^] E ^ [yy 



Tl 



E [yx"^] < E [xx-^] . Note that the inequality in (jl44p is strict because of the linear independence of w^'s-D 



I Proof of Proposition 138 



Statement (a) follows from Proposition I30t for (b) and (c), we prove only (b). 

According to Proposition 128(a) . we need to show that, under the assumption that T™° in ()74p is positive 
semidefinite, the augmented Jacobian matrix JF'^(Q) associated to F*^(Q) in (j73p is augmented positive 
semidefinite on , Given ([73|) . Dq*F'^(Q) = 0, implying that JF'^(Q) is a block diagonal matrix: 

Z)qFC(Q) 

(I?qFC(Q))* 

with Z)qF''-'(Q) given by (see Appendix [Fl for a similar computation): 

i^Q,Ff(Q,) ••• DQ.FfiQ) 



JF^(Q) = - 



^qF^(Q) 



Z)q,F^(Q) ... Da^FfiQ) 

*ii(Q)*®*ii(Q) ••• *i/(Q) 



*n(Q)*»*n(Q) 



*i/(Q) 
*//(Q)*®*//(Q) 



with 



*jj(Q) — Hf R„. + ^HjjQjHj^ Hij — H^Si(Q)Hij. 



(145) 



(146) 



(147) 



(148) 



=S.(Q) 
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1 /2 

We will denote by (Q) the square root of the positive definite matrix ^'jj(Q) (recall that the channel 

U /2 1 /2 

matrices tin are assumed to be full-column rank), i.e., ^'jj(Q) = (Q) (Q). 

Therefore, JF*^(Q) is augmented positive semidefinite on 'P'^^'^° if and only if I?qF^(Q) is positive semidef- 
inite on 'p^^'^° ^ or equivalently the following matrix is so: 



1/2 



H/2 



17*11 



-1/2 



^11 ^n^H 



n^ii 



-1/2 



(149) 



where for notational simplicity we omitted the dependence on Q and write ^ij, instead of ^'jj(Q). The 
condensed matrix associated to (1149^ is the following I x I matrix 



-1/2 



II ^Il^H 



-H/2 



II 



-1/2 



(150) 

where || A||2 — \/ p{A^ A) is the spectral norm of A. Note that we can rewrite each of the off-diagonal terms 
of (dSOD as: with * 



-1/2 



^11 ^ij^ii 



T~* ~7f~ mV2 



:i5i) 



where in the last equality we used the property p (A-^A* (8) A^A) = p (A'^A*) p (A^A) and the fact that 
the eigenvalues of A-^A* coincide with those of A-'^A. Using (jlSip . we can now introduce the so-called 
comparison matrix T™°(Q) associated to (|150p and defined as 



pi^M) *.,(Q) 



if i = j 
otherwise. 



It is indeed not difficult to see that if Tpi"^°(Q) is positive semidefinite on 



then so is the matrix (jl50p and 
thus also DqF'^(Q). To conclude the proof, it is enough to show that T^^'"°(Q) > "^fc" for all Q G -p^i""", 
where Tp^(?° is defined in (I74p and the inequality has to be intended component-wise. The latter properties 
indeed implies that if T™° is positive semidefinite then so is T™°(Q) on 'P"'^'"°. To this end, we focus 



^ I ['Y^mimoi I 



next on the off-diagonal terms /)(*ij- (Q) *ij(Q)) of T^c"°(Q) and prove that |[T^"c"°(Q)]iil 

1 /2 

for ah Q G P""^""" and i / j. Denoting by St'"(Q) the square root of the positive definite matrix Si(Q) 



defined in (fTigl) [i.e., Si(Q) 



(Q)S,'/'(Q)], and using *,,(Q) ^ *r-/^(Q)*,^.(Q)*-^/^(Q), we have 



HI2, 



-1/2/ 
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the following chain of equalities/inequalities: for all Q E "puinio and i ^ j, 

= p (*,::^/'(Q)H,(^S,(Q)H,, (Hf S,(Q)H,,)"'Hf S,(Q)H,,*ri/2(Q 
= p (*,::^/^Q)H,f Sf/^Q) Sp(Q)H,, (Hf S,(Q)H,,)"' Hf sf^^Q) Sp(Q)H,,*,::^/'(Q) 



(152) 

< p (sj/^Q)H,,*,::nQ)Hf Sf /'(Q)) = p (sj/^Q)H,, (Hf S.(Q)H,,)-' sP(Q)) 

(153) 

< p (s,^/^(Q)H,,Ht^Sri(Q)Htf Hf^sP(Q)) (154) 



< p (S-i(Q)) • p(S,(Q)) . p (njf H^H,,hO (155) 



<INNR,,.p(HlfHfH,,Hl, 



[T^r],, (156) 



where in (fT^ . P7^(h,,) = Sj/^(Q)Hii (HfSi(Q)Hii) ^Hfsf/^(Q) is the orthogonal projection onto the 
range space of Hjj; in (jl53p we used the property of the projection P7^(Hi,) ^ I ^"^^ the spectral radius 
inequality p(A^BA) < p(A^CA) for all ^ B ^ C; (fTSD foUows from the property (X^AX)"^ < 
X^A^^X^''^, valid for all positive definite n x n matrices A and n x k full-column rank matrices X; in 
(|155p we used A ^ P ( A) • I and the spectral radius inequality as in (|153|) : and finally (|156p follows from 
p(Sri(Q)) •p(Si(Q)) < INNRjj, with INNR^j defined in ([751) . 

The above chain of inequalities proves the desired relationship between Tp^j?°(Q) and T™°, which 
completes the proof. □ 
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